Functional Comparison of Innate Immune Signaling Pathways in Primates

Summary

Humans respond differently than other primates to a large number of infections. Differences in susceptibility to infectious agents between humans and other primates are probably due to inter-species differences in immune response to infection. Consistent with that notion, genes involved in immunity-related processes are strongly enriched among recent targets of positive selection in primates, suggesting that immune responses evolve rapidly, yet providing only indirect evidence for possible inter-species functional differences. To directly compare immune responses among primates, we stimulated primary monocytes from humans, chimpanzees, and rhesus macaques with lipopolysaccharide (LPS) and studied the ensuing time-course regulatory responses. We find that, while the universal Toll-like receptor response is mostly conserved across primates, the regulatory response associated with viral infections is often lineage-specific, probably reflecting rapid host–virus mutual adaptation cycles. Additionally, human-specific immune responses are enriched for genes involved in apoptosis, as well as for genes associated with cancer and with susceptibility to infectious diseases or immune-related disorders. Finally, we find that chimpanzee-specific immune signaling pathways are enriched for HIV–interacting genes. Put together, our observations lend strong support to the notion that lineage-specific immune responses may help explain known inter-species differences in susceptibility to infectious diseases.

Introduction

Due to our natural focus on humans, we know of a large number of diseases or medical conditions that affect humans more severely than non-human primates. Examples include progression to AIDS following infection with HIV, progression to malaria following infection with Plasmodium falciparum, Alzheimer's disease, cancer, and adverse complications following infection with hepatitis B and C (reviewed in [1], [2]). Differences in susceptibility to infectious agents between humans and other primates might be explained, at least in part, by inter-species differences in immune response to infection. Indeed, a large body of work indicates that immune systems are rapidly evolving. In particular, while very little comparative functional data in primates has been collected, recent genomic scans for signatures of natural selection have reported that genes involved in immunity processes are strongly enriched among targets of positive section in human and chimpanzee [3]–[11].

Immune responses are typically classified as either ‘innate’ or ‘adaptive.’ Historically, the focus of most immunological studies has been on the adaptive response and its hallmarks, namely the generation of a large repertoire of antigen-recognizing receptors and immunological memory. Recently, however, more effort has been expended on understanding the innate immune system, as it became clear that innate immunity is an evolutionarily ancient defense mechanism, which governs the initial detection of pathogens and stimulates the first line of host defense [12]–[15]. Moreover, innate immune responses were shown to play a pivotal role in the development of pathogen-specific humoral and cellular adaptive immune responses, which are mediated by B and T cells [16]–[18].

The recognition of pathogens by the innate immune system is primarily mediated by phagocytic cells (e.g., monocytes, macrophages, and dendritic cells) through germline-encoded receptors, known as pattern recognition receptors (PRRs) [17], [19], [20]. The PRRs recognize conserved molecular features characteristic of the microbial world, commonly referred to as pathogen-associated molecular patterns (PAMPs) [17], [19], [20]. Among the different PRRs, the Toll-like receptor (TLR) family, which comprises 10 functional members in humans, has been the most extensively studied. For example, by stimulating primary cell cultures with different TLR agonists in vitro (e.g., references [21]–[23]) and by studying mouse models that lack one or several TLRs (e.g., references [24]–[26]), it has been shown that TLRs can be activated in response to virtually any microbe that invades the host.

Once activated, TLRs play a crucial role in orchestrating the response to pathogenic microbial infections through the induction of two major regulatory programs. First, a universal regulatory response, which can be activated by all TLRs and is triggered by infection with a diverse range of microbes or TLR agonists [21], [22], [27]–[29]. This response has been interpreted as a generic ‘alarm signal’ for infection [22], [28], [29]. Second, individual TLRs can activate regulatory programs that are specific to individual microbial agents [19]. Comparative functional studies of TLR-mediated immune response in primates might therefore shed light on inter-species differences in susceptibility to certain infectious agents. However, at present, there is very little functional data with which one can study the evolution of the immune system in primates.

Results

In order to study functional differences between the innate immune response of humans and two close evolutionary relatives, chimpanzees (Pan troglodytes) and rhesus macaques (Macaca mulatta), we stimulated primary monocytes from six individuals from each of the three species with LPS for 4, 12, and 24 hours (see Figure S1 for an illustration of the study design). LPS activates the TLR pathway (specifically, TLR4) and mimics an infection with Gram negative bacteria [16]–[18]. We chose this treatment because LPS, via TLR4, activates multiple immune signaling pathways, leading to the induction of both inflammatory and ‘viral-like’ responses [19]. Additionally, the stimulation of immune cells with LPS was shown to result in a very similar regulatory response (87% overlap) to the response to infection with a live bacteria such as E. coli [22].

To confirm that the LPS treatment activated TLR4-mediated immune responses, we used quantitative PCR to estimate the induction levels of three inflammatory cytokines (IL6, IL1β, and TNF). In all samples (from all individuals at all time points), levels of the three inflammatory cytokines were significantly higher following stimulation with LPS (Figure S2). However, we noticed that the quantitative responses to the treatment in the chimpanzee samples were lower than those of the human and rhesus macaque samples. This observation probably reflects a technical difficulty in culturing chimpanzee primary monocytes without inducing a general stress response, which results in the attenuation of the quantitative response to further stimuli (see Materials and Methods for more details). In what follows, we therefore focus primarily on qualitative rather than quantitative differences between individuals and species in the regulatory response to stimulation with LPS.

LPS-mediated immune responses in primates

To estimate and compare gene expression levels in samples from multiple species, we used a multispecies microarray, which includes orthologous probes from human, chimpanzee, and rhesus macaque for 18,109 genes [30]. Following processing and normalization of the array data, we used a gene-specific linear mixed-effect model (see Materials and Methods) to identify inter-species differences in the regulatory response to stimulation with LPS (the ‘treatment’). To minimize the number of falsely identified differences across species, we applied two statistical cutoffs for classifying genes as responding to the treatment. Specifically, conditional on observing a treatment effect with high statistical confidence in one species, we assumed that a treatment effect likely occurred in other species as well, and relaxed the statistical cutoff for the classification of such secondary observations (see Materials and Methods for more details and the specific statistical cutoffs used). This procedure minimizes the number of falsely identified inter-species differences that might ultimately arise from incomplete power to identify differences in gene expression levels following the treatment.

Using this approach, we identified 3,170 genes whose expression levels changed following the treatment in at least one species, at any time point, of which 793 genes responded in all three species (Figure 1A, Table S1, Figure S3). As expected, genes that responded to stimulation with LPS in all three species are enriched with genes involved in immune-related biological processes such as “inflammatory-response” and “cytokine-signaling”, as well as in specific immune-related pathways including the Toll-like receptor pathway, cytokine-cytokine receptor interactions, and the Jak-STAT signaling pathway (FDR for all reported results is <0.01; Figure 1B, 1C, Table S2). Consistent with previous observations in functional studies of the immune system in mice [21], [27], we found that the conserved regulatory response to stimulation with LPS in primates included an enrichment of genes that are likely regulated by the transcription factor NF-kB (P<10⁻⁵) and several interferon regulatory factors (e.g., IRF7 and IRF1; P<10⁻³; Figure 1D, Table S3). Put together, these observations clearly demonstrate that the monocytes from all three species responded to the treatment with LPS by engaging TLR4-mediated regulatory pathways [19], leading to the induction of pro-inflammatory and anti-viral immune responses via the activation of NF-kB and IRF mediated pathways.

To gain further insight into the evolution of LPS-induced immune responses in primates, we classified genes as participating in either the universal regulatory response to infection (which can be triggered by a diverse range of microbes or TLR stimuli [21], [22], [27]), or in the microbial-specific response (which we then further classified as responses to either bacterial or viral infections [21], [22], [27]). Based on these classifications (Table S4), we examined how genes falling into each of the categories responded to infection in humans, chimpanzees, and rhesus macaques. We found that the majority (58%) of genes involved in universal response to infection showed a conserved regulatory response to stimulation with LPS in all three species, compared to only 31% of genes known to respond primarily to either viral or bacterial infection (χ² test, P<0.001; Figure 2A). Viewed from a different perspective, we observed that the proportion of genes involved in immune response to viral infections is significantly higher (1.5-fold) among genes that responded to stimulation with LPS in only one of the species, compared with genes that responded to the treatment in all three species (χ² test, P = 0.002; Figure 2B). Taken together, our data strongly support the notion that the universal TLR response is mostly conserved across primates and that much of the divergence in immune response is observed in genes that are involved in response to specific microbial and viral agents.

Characterization of species-specific innate immune responses

We proceeded by focusing on species-specific immune responses. Using the conservative approach described above, we identified 335, 273, and 393 genes as responding to stimulation with LPS exclusively in the human, chimpanzee, and rhesus macaque monocytes, respectively (see Figure 3A–3C for examples). To characterize these gene sets, we considered functional annotations based on the GO and KEGG databases (Table S5, S6, S7). Somewhat surprisingly, the only significant enrichments (after correction for multiple tests) were observed among the 335 genes that responded to the treatment exclusively in humans. We found that human-specific immune response was enriched for genes in pathways previously associated with cancer (e.g., Chronic myeloid leukemia or prostate cancer; P≤3.0×10⁻³, FDR<0.06), the B cell receptor signaling pathway (P = 3.2×10⁻³, FDR = 0.06), and pathways related to apoptosis (P = 5.0×10⁻³, FDR = 0.07; see Table S5 for a complete list of significant results). Further, by using the STRING database [31] to visualize all known functional interactions between these 335 genes, we found that 151 of the genes in this set (45%) are known to interact with each other – using the default cutoff suggested by STRING to define a functional interaction (Table S8). Applying a more stringent cutoff (a STRING confidence-score higher than 0.7), we identified 78 genes (23%) that interact with each other, in a functional module that is enriched with genes involved in cancer biology and apoptosis pathways (Figure 4, Table S9). In order to obtain further support for interactions across these 78 genes, we used GRAIL, a tool that uses text mining of PubMed abstracts to identify published functional interactions between genes. We found that 43 out of the 78 human-specific immune response genes (55%) had a GRAIL score of P-text<0.05, statistically supporting the notion that they have a functional interaction with at least one other gene in the list (only ∼7% of genes are expected to have GRAIL score of P-text<0.05 in randomly chosen sets of 78 genes).

We then considered networks of co-expressed genes (namely genes with coordinated patterns of expression) for each species, to find additional putative modules of interacting genes (see Materials and Methods). We found 33, 17 and 32 regulatory modules in humans, chimpanzees and rhesus macaques, respectively, with an average connectivity (|r|) higher than 0.5 (Figure 5, Table S10, S11, S12). Based on 100 random permutations of the gene expression values, we estimated that the number of clusters with |r|>0.5 expected by chance alone is 1.28±1.04, 1.16±1.08, or 1.37±1.08, using data from humans, chimpanzees and rhesus macaques, respectively, suggesting that the observed excess of regulatory modules likely describe meaningful biological relationships.

In humans, the largest regulatory module contains 82 genes, which are significantly enriched for several biological processes involved in extracellular matrix remodeling (P<0.002; Table S10). The second largest regulatory module is significantly enriched for genes involved in apoptotic pathways (P<0.04; Table S10). Interestingly, “apoptosis-related” processes appeared to also be enriched among genes in one of the largest regulatory modules identified in chimpanzees, with 23 genes (Table S11), as well as among three large regulatory modules (>30 genes each) identified in rhesus macaques (Table S12). As the regulatory modules comprise of mutually exclusive sets of genes across the three species (by the nature of the analysis), these observations support the notion that immunological-associated apoptosis mechanisms evolve rapidly in primates.

Finally, we asked whether the observed inter-species differences in immune response might provide insight into the mechanisms underlying differences in susceptibility to infectious diseases between humans and non-human primates. To do so, we considered the subsets of genes that responded to stimulation with LPS exclusively in humans, chimpanzees or rhesus macaques, and examined whether they were enriched for genes previously reported to be associated with immune disorders in humans and/or susceptibility to infectious diseases (see Materials and Methods).

We found an enrichment of “immune-related-disease-genes” among genes that responded to the treatment with LPS exclusively in humans (χ² test, P = 0.03; Figure 6A). Interestingly, we also found that the set of genes that responded to stimulation with LPS exclusively in chimpanzees was enriched with genes that code for host cell proteins known to interact with HIV-1 (Figure 6B; χ² test, P = 0.0002). No significant enrichment of HIV-1 interacting genes was observed among genes that responded to stimulation with LPS exclusively in either humans or rhesus macaques. This observation is robust with respect to the specific cutoffs used to classify genes that responded to stimulation in LPS in only one species (Figure S5).

Discussion

We have performed a genome-wide study of LPS-mediated immune responses in primary monocytes from humans, chimpanzees, and rhesus macaques. Our study design allowed us to characterize conserved innate immune response mechanisms in primates as well as to identify species-specific regulatory responses to stimulation with LPS.

An important difficulty of all studies of gene regulation in primary tissues from primates, apes in particular, is the inability to stage the environment for each of the donor individuals across species. In our study, biological replication within species partially addresses this difficulty, but the possibility that a subset of the observed inter-species differences in gene regulation are due to differences in environments (e.g., diet) across species still exists. An additional difficulty is that most available tools for manipulating cell cultures and performing immune-related assays have not been optimized to work with non-human primate cells. We addressed this issue by performing a large number of quality controls, including the validation of the response to stimulation with LPS in each cell culture by using qPCR. Nevertheless, our observation of systematic differences in the quantitative response to infection across cultures from different species (in particular, from chimpanzees) probably has a technical rather than a biological explanation. For that reason, we chose to draw conclusions primarily based on qualitative differences between species. Thus, the inter-species regulatory differences reported in our study likely provide a lower boundary for the actual number of differences in immune response between humans, chimpanzees and rhesus macaques.

Conserved innate-immune responses in primates

We identified 793 genes that responded to stimulation with LPS in all three species. As expected, this set of genes was significantly enriched for genes involved in immune responses, and specifically for genes involved in TLR-mediated pathways. Some examples of conserved TLR4-induced immune responses include the strong up-regulation of several pro-inflammatory cytokines, such as IL-6, IL1-β and tumor necrosis factor (TNF), and chemokines, such as CCL2, CCL3 and CCL4, whose roles are to recruit other effector cells to the inflammatory site [29], [32]. We also observed a conserved up-regulation of the anti-inflammatory cytokine IL-10, probably to control the levels of inflammatory response and avoid tissue damage [29], as well as the up-regulation of several interferon-α inducible genes (e.g., IFIH1, IFIT1, and IFIT3).

Overall, conserved immune responses were enriched for genes whose expression levels are regulated by the transcription factor NF-kB, or by several interferon regulatory factors, which are the master regulators of TLR4-dependent pathways. Interestingly, before infection, the expression levels of many of these master regulators (e.g., REL, NFKB1, RELB, IRF2, IRF9) were different across the three species, while post-infection, their expression converged to practically the same level, regardless of species (Figure S6). This observation suggests that the regulatory response of these key transcription factors likely evolve under strong evolutionary constraints, probably to ensure efficient downstream immune responses.

A known property of the regulatory programs mediated by different TLRs is the activation of both a universal response (shared by all TLRs) as well as a response that is specific to each microbial agent (or TLR ligand) [21], [22], [27]–[29]. We found that the universal TLR response is remarkably more conserved across primates compared to microbial-specific responses. From an evolutionary perspective this observation makes intuitive sense. Indeed, ‘core’ immune responses, which are critical to fight any invading pathogen, are expected to be under stronger evolutionary constraint compared to immune programs that are only important in the presence of specific microbial infections. Consistent with this expectation, our data also support the notion that adaptation of innate immune responses in primates primarily took place at the level of ‘peripheral’ responses, namely, pathogen-specific immune responses.

Species-specific immune responses

Among genes whose regulation was affected by stimulation with LPS in only one species, we found an enrichment of genes associated with response to viral infections. This observation might reflect the need of the host immune system to frequently devise new defense mechanisms to fight viral infection, as viruses tend to evolve faster than other microbes [33]. We also found that species-specific immune responses are enriched with genes annotated to have a role in apoptotic pathways. Apoptosis is a critical component of successful immune response as infected cells have to be efficiently removed without inciting an inflammatory reaction [34]. Moreover, controlled cell death is used to restore normal cell numbers following clonal expansion of antigen-specific lymphocytes [34]. Consistent with our observation of rapid evolution of the regulation of apoptotic pathways, coding regions of apoptosis-related genes have previously been shown to be rapidly evolving during primate evolution [7], [35]. Put together, these observations suggest that inter-species differences in apoptosis-related immune responses may be adaptive. Although the selective pressures underlying these adaptations are unclear, these observations might help elucidate the basis for important phenotypic differences between humans and non-human primates, such as differences in susceptibility to cancer.

Indeed, cancer incidence in non-human primates is low compared to that observed in humans, even when age is taken into account [1], [2], [36]–[39]. The deregulation of apoptosis has been extensively described as a hallmark of cancers [40]. Thus, while our observation that the human-specific immune response to stimulation with LPS is characterized by a significant enrichment of cancer-related genes is not surprising (because the ‘cancer-related’ and ‘apoptosis-related’ gene sets are not mutually exclusive), it may provide a first step towards understanding the mechanisms underlying the differences in cancer incidence between humans and other primates. For example, we observed that the pro-apoptotic gene CASP10 was strongly down regulated early after stimulation with LPS, exclusively in humans (Figure S7). Somatic mutations in CASP10, as well as reduced expression levels of this gene, were found to be associated with a number of different human cancers [41]–[44]. The observed inter-species differences in the regulation of CASP10 following infection may therefore be related to differences in the rates of cancer across species. Detailed comparative studies of apoptosis-related regulatory mechanisms in model organisms will be necessary to fully explore the possible connection between the predisposition to cancer and inter-species differences in immune responses.

Inter-species differences in susceptibility to HIV/AIDS

Genes whose regulation was altered following stimulation with LPS exclusively in humans were enriched with genes known to be associated with susceptibility to infectious diseases or to immune-related diseases. We did not observe such enrichment when we considered the immune responses specific to chimpanzees or rhesus macaques. Our observations make intuitive sense, as we know more about the genes associated with diseases that affect humans than those that affect the two non-human primate species. In other words, it is reasonable to assume that we would have found similar enrichments in chimpanzees and rhesus macaques if we knew more about the genetic basis of infectious and immune-related diseases that primarily affect these two species. Our observations thus underscore the link between species-specific immune responses and susceptibility to infectious disease.

One interesting example is the enrichment of genes known to interact with HIV among genes whose regulation was affected by the stimulation with LPS exclusively in chimpanzee. This observation is intriguing because, unlike humans and rhesus macaques, a large number of studies propose that chimpanzees only rarely develop AIDS following infection with HIV [45]–[47]. This notion has recently been challenged by Keele et al. [48], who reported that wild chimpanzees naturally infected with SIVcpz do develop hallmarks of AIDS. The apparently contradictory observations in the literature might be explained by the fact that Keele et al. used data collected from an eastern subspecies of chimpanzees (Pan troglodytes schweinfurthii), whereas previous observations of increased protection from AIDS, were based on studies with the western chimpanzee subspecies (Pan troglodytes verus), the one used in our study. Differences in susceptibility to HIV between sub-species of chimpanzees might be explained by the fact that Pan troglodytes verus were physically separated from the two other sub-species prior to systemic infection with the two recombinant monkey viruses of SIVcpz.

Some examples of HIV-interacting gene that responded to stimulation with LPS exclusively in chimpanzees include ITGB2(CD18) and ITGAM(CD11b), which are the two members of the complement receptor 3 (CR3) that have been shown to play a key role in the infection of dendritic cells by C3-opsonized HIV [49], [50] and the viral transfer to CD4 T cells [50]. Interestingly, these two genes were down-regulated after LPS stimulation, exclusively in chimpanzee monocytes. Another example is the APOBEC3F gene, which is one of the most potent inhibitors of HIV replication [51] and was significantly up-regulated in response to stimulation with LPS, only in chimpanzee monocytes. The direction of regulatory change in these cases (namely, the down regulation of a receptor that may be used by the HIV virus, and the up-regulation of a known inhibitor of HIV replication), is consistent with a theoretical mechanism of increased resistance of chimpanzees (or at least Pan troglodytes verus) to progression of AIDS. That said, future studies are now required to evaluate if the down-regulation of CR3 and up-regulation of APOBEC3F are also observed after infection with HIV (or SIVcpz).

Our observations may reflect an adaptation of the chimpanzee immune system to infection with HIV/SIV or perhaps to other retroviral infection(s). Previous studies of variation at the nucleotide level have reported that genes associated with HIV infection (such as CD45, APOBEC3G and APOBEC3H) evolved under positive selection in primates [52], particularly after the divergence of humans and chimpanzees [11]. Taken together, our data suggest that regulatory changes occurring specifically in the chimpanzee lineage might explain, at least in part, why chimpanzees tend not to progress to AIDS following infection with HIV/SIV.

More generally, our observations may help to explain other inter-species differences in susceptibility to infectious agents, such as the increased resistance of chimpanzees to certain other viral infections, including hepatitis B and C, and influenza A. Our study, however, is only the first step in characterizing inter-species differences in immune response, in particular because LPS is a general stimulant. We expect future comparative studies in primates to focus on the immune response to different individual infectious agents.

Materials and Methods

Ethics statement

This study was conducted according to the principles expressed in the Declaration of Helsinki. An IRB approved consent form was obtained from each human donor. Collection of the non-human primate samples was perform at the Yerkes Primate Center, in a manner that conformed to the animal subject regulatory standards enforced by the Emory University Institutional Animal Care and Use Committee (IACUC approved protocol #028-2009Y).

Sample collection

We measured gene expression levels in blood monocytes from six humans, six chimpanzees and six rhesus macaques (three males and three females from each species; see Table S13 for details on all samples). Blood samples were collected in BD Vacutainer CPT Cell Preparation Tube (BD, Franklin Lakes, NJ) and peripheral blood mononuclear cells (PBMCs) were purified according to the manufacturer's instructions. Non-human primate blood samples were collected at the Yerkes primate center and human samples were obtained from Research Blood Components.

Monocyte purification and culturing

Blood monocytes from the three species were purified from PBMCs using magnetic cell sorting technology (MACS technology from Miltenyi Biotech). Specifically, monocytes from humans and rhesus macaques were purified by positive selection with magnetic CD14 MicroBeads (Miltenyi Biotech). This method did not work well with the chimpanzee samples (less than 3% of chimpanzee PBMCs were isolated using cell sorting with a CD14 antibody). Instead, monocytes from chimpanzees were purified by depletion of non-monocyte cell types using the “Monocyte Isolation Kit II” (Miltenyi Biotech). Regardless of the method used, the purity of the isolated monocyte population was evaluated by flow cytometry. For the human and rhesus macaque samples we used a fluorochrome-conjugated antibody against monocytes (CD14-FITC; Beckmancoulter). For chimpanzees, we further confirmed the purity of the monocyte population by also using two additional fluorochrome-conjugated antibodies - against B-cells (CD20-PE; BD Bioscience) and T-Cells (CD3-APC BD Bioscience). Regardless of species, only samples with monocyte purity higher than 80% were used in subsequent experiments (Figures S8 and S9). Finally, to further increase the purity of the monocyte-enriched fraction (virtually to 100%), we performed an additional selection step based on the unique capacity of monocytes/macrophages to strongly adhere to the plastic of cell culture dishes. To do so, we cultured the cells overnight (see below for details) and washed the cell culture wells in the morning to retain only adherent cells (i.e., monocytes).

Since we used different methods to purify monocytes in samples from chimpanzees (negative selection) than in samples from humans or rhesus macaques (positive selection), we also performed a control experiment to empirically evaluate to what extent the purification method used affects the ability of the purified cells to respond to LPS stimulation. To do so, we collected whole blood samples from three additional humans, from which we purified monocytes using both a positive and a negative selection method. We then performed the same LPS treatment experiment we applied to the main samples from all three species (described below), and compared differences in response to LPS treatment between monocytes purified by positive selection and monocytes purified by negative selection. As discussed in Text S1 and Figures S15 and S16, we found that the method used to purify the monocytes (negative or positive selection) has only a minimal effect on the measured regulatory response to stimulation with LPS. That said, to be conservative, we excluded from all analyses presented in the manuscript the 192 genes identified as responding to LPS treatment only in monocytes purified by either positive or negative selection.

Monocyte culture and LPS stimulation

Monocytes were cultured in 24-well cell culture plates (Corning) in serum free media (CTL's test media) at a density of 1 million cells per ml. We used a serum free media to minimize the probability that the monocytes were non-specifically activated as a result of the undefined nature of serum products (e.g., as a result of a mitogenic serum batch). The cells were then stimulated with 1 ug/ml of LPS (Invivogen, Ultrapure LPS, E. coli 0111:B4) for 4, 12, and 24 hours. These time points were chosen based on previous observations that the transcription kinetics of immune response to infection can generally be characterized by early, middle, and late phases of response, which can be effectively captured at 4, 12 and 24 hours post infection [22]. All time course experiments were started at the same time of day (∼8am) to prevent the introduction of variation due to differences in circadian rhythm. Untreated cell cultures were kept alongside the stimulated cultures and were harvested at the same time intervals (4, 12, and 24 hours post stimulation).

RNA extraction and amplification

Total RNA from each cell culture was extracted using RNeasy columns (Qiagen, Valencia, CA). For all samples, RNA quantity was evaluated spectrophotometrically, and the quality was assessed with the Agilent 2100 bioanalyzer (Agilent Technologies Inc, Palo Alto, CA). Only samples with no evidence for RNA degradation (RNA integrity number >8.5) were retained for further experiments. To evaluate the activation of monocytes after stimulation, we used quantitative PCR to test for an over-expression of Tumor-necrosis Factor (TNF-a), Interleukine-6 (IL-6) and IL1-β. inflammatory cytokines that are known to be induced following the activation of TLRs (primer sequences and PCR conditions can be found in Table S14). Only samples for which we observed a significant induction of these cytokines were used in downstream experiments. Once we confirmed that monocytes from all three species responded to the treatment, we performed linear amplifications of the total RNA samples by using in-vitro transcription. Specifically, 400 ng of high-quality total RNA were amplified using the MessageAmp II kit (Ambion). Unlike exponential RNA amplification methods, aRNA amplification has been shown to maintain the relative representation of the starting mRNA population [53], [54] (Figure S10).

Multi-species microarray

To compare genome-wide gene expression levels between humans, chimpanzees, and rhesus macaques, we hybridized the RNA samples to the multi-species microarray described by Blekhman et al. [30]. This array contains orthologous probes from the three species, thus allowing a comparison of gene expression levels between species without the confounding effects of sequence mismatches on hybridization intensities [30]. The microarray contains probes for 18,109 genes (see Blekhman et al. [30] for a detailed description of the multi-species array). The labeling of the amplified RNA samples and subsequent hybridization to the microarray were performed by Nimblegen. For each individual we hybridized one non-stimulated and one stimulated sample at each of the three time points (4 hours, 12 hours and 24 hours). The total number of arrays analyzed was therefore 108 ( = 3 species ×6 individuals ×6 arrays per individual). Quality control, background correction and normalization of the expression data were performed as previously described [30] (Figure S11 and S12).

Statistical analysis

All the statistical analyses detailed in this and the following sections were performed using the R statistical environment (http://www.r-project.org).

Identifying genes differentially expressed between species

We analyzed data from the 17,231 genes (95% of the genes on the array) that were assayed by at least three orthologous probes across all species. To identify differentially expressed genes between the three species, across time-points and following the treatment, we modeled the expression levels of each gene independently by using a linear mixed-effects model similar to the one described in Blekhman et al. [30]. Specifically, for each of the 17,231 genes, if y_sroi denotes the normalized log₂ intensity expression value for individual i (i = 1,..6), from species s (s = human, chimpanzee or rhesus macaque), measured at probe r (r = 1, …,7), which is derived from species o, we assume that:(1)where:

Here, µ_s is a species-specific fixed-effect (representing changes in expression levels across the three species), π_ro is a fixed-effect representing the probe effect for each individual probe within a probe-set and the composition effect of species-specific orthologous probes, and κ_sro is a fixed-effect representing the attenuation of hybridization intensities due to sequence mismatches between species of RNA and a species-specific derived probe, which are different for each individual probe within a probe-set (see reference [30] for more details). The term represents a random-effect for individual i from species s; this effect is assumed to follow a Normal distribution with mean 0 and variance . To determine whether a gene was differentially expressed between species, we assessed how well model (1) fitted the data under the following parameterizations of :

In the null model () the gene's expression level is assumed to be constant across the three species, while in the alternative model () the expression level is allowed to differ between species. We maximized the likelihood under these two parameterizations using an iterated least-squares approach and compared the fit of the models by calculating the likelihood-ratio test statistic. We calculated p-values based on a χ² distribution with two degrees of freedom, and corrected for multiple testing using the FDR approach of Benjamini and Hochberg[55].

Identifying genes whose expression level was altered following stimulation with LPS treatment

Due to the paired nature of the study design (untreated and treated cells from the same individual are compared), parameterizations of an extended version of model (1) are difficult to interpret. Instead, we first regressed out probe effects and then used a linear model framework to identify genes whose expression levels have changed following the treatment. To regress out probe effects, we considered data from each time point separately. For each gene, if y_sroit denotes the normalized log₂ intensity expression value for individual i (i = 1,..6), from species s (s = human, chimpanzee or rhesus macaque), from class t (t = non-stimulated or stimulated) measured at probe r (r = 1, …,7), which is derived from species o, we assume that:(2)where:

In the above, , , and are defined as in the previous section. However, we have added two additional terms to the model: (a treatment effect) and (a treatment-by-species interaction). Further, the parameterization of the random effect has been modified to reflect the incorporation of the treatment effect, such that . We used a maximum likelihood approach to estimate the parameters in this model, and subsequently calculated a corrected measure of expression, for each individual i from each species s and for the two classes t as: ( denotes the maximum likelihood estimate of et cetera). We analyzed these corrected measures of expression levels further using gene-wise models. Specifically, for each species separately, if denotes the corrected gene expression level in individual i undergoing treatment t at the four hour time-point, we fitted the following two models:

Here, corresponds to an intercept, corresponds to an individual effect (shared across treatments), corresponds to a treatment effect and is assumed to follow a normal distribution with mean 0 and variance . We tested whether there was significant evidence of a difference in expression levels between the treated and un-treated samples by comparing the fit of these models using a likelihood-ratio test statistic. We calculated p-values based on a χ² distribution with one degree of freedom.

Ultimately, our goal was to identify inter-species differences in the regulatory response to stimulation with LPS. To minimize the number of falsely identified differences across species we applied two statistical cutoffs for classifying genes as responding to the treatment. Specifically, as a first step, at each time point, using data from the human and rhesus macaque samples, we initially classified genes as differently expressed following the treatment using an FDR cutoff <0.001, also requiring that the effect size of the treatment was equal or larger than |0.3| (log₂ scale). For data from the chimpanzee samples, we chose to use different initial cutoffs (an FDR<0.1, and an absolute effect size cutoff of 0.1), as the effect sizes for the regulatory response to the treatment were overall lower (importantly, however, the overlap in lists of responding genes across species is high regardless of the specific cutoff used; Figure S13, Table S15). We then assumed that, conditional on observing a treatment effect with high statistical confidence in one species, a treatment effect likely occurred in other species as well. To classify such secondary observations we used a relaxed cutoff of P<0.05. This approach is therefore conservative with respect to identifying differences in immune response across species. While the specific numbers reported in the paper are based on arbitrarily chosen cutoffs, the qualitative conclusions we discuss are robust with respect to the specific choice of cutoffs within a considerable range (Table S16 and Figure S5).

Gene ontology and pathway enrichment analysis

We used GeneTrail (http://genetrail.bioinf.uni-sb.de) [56] to test for enrichment of functional annotations among different classes of genes (as detailed in the results). In all tests, we used a background set of 13,244 genes, which were classified as expressed (using an absolute log₂ intensity cutoff of 7.5; Figure S14) in at least one species at one condition; that is, any of the time points for the non-stimulated or stimulated samples. The tests were performed using all GO categories and KEGG pathways. We calculated p-values using a Hyper-geometric distribution, and used the approach of Benjamini and Hochberg [55] to control the false discovery rate.

Enrichment of transcription factor binding sites

We applied the promoter analysis algorithm PRIMA implemented in the EXPANDER package [57]. Given a target set and a background set of genes, PRIMA identifies transcription factor binding motifs that are significantly more prevalent in the promoter of the target set than in the background set. As background we used all the genes in the array that had at least 3 homologous probes across species and that were expressed in at least one species in one of the six conditions (i.e., any of the three time-point in either treated or non-treated samples).

Network analysis

We looked for known functional associations between the 335 human-specific LPS response genes using the STRING database (http://string.embl.de/). STRING is a database of both known and predicted protein-protein interactions, which includes direct (physical) and indirect (functional) associations derived from numerous sources, including experimental repositories, computational prediction methods and public text collections [31]. We selected all interactions/associations available for a given node with a combined score greater than 0.7. The scores given in the STRING database define the confidence limit for each described interaction/association. A combined score of 0.7 is recommended as the high stringency criterion by the database authors.

In turn, to identify modules of co-expressed genes, we used the Modulated Modularity Clustering (MMC) algorithm, which seeks community structure in graphical data; that is, a graph of genes connected by edges whose weights reflect the degree to which their transcriptional profiles are correlated (see reference [58] for a detailed description of the method). Co-expression modules were defined using the probe-corrected expression estimates for the set of 335, 273 and 393 genes whose expression levels were altered following the treatment exclusively in humans, chimpanzees and rhesus macaques, respectively.

Analysis of genes associated with disease

Genes previously reported to be associated with immune disorders and/or susceptibility to any infectious disease were identified using the Genetic Association Database (http://geneticassociationdb.nih.gov/). We downloaded the full GAD dataset on Nov 9, 2009, and parsed the all.xls table, which contains gene-disease associations. We labeled genes as immune-related if they had been associated with diseases for which the ‘disease class’ field in the GAD database was defined as ‘Immune’ or “Infectious diseases”. The list of host genes known to interact with HIV-1 proteins was retrieved from the HIV-1 Human Protein Interaction Database, which catalogues over 1,400 human proteins reported in the scientific literature to participate in HIV-1 to human protein interactions [59].

Supporting Information

Attachment 1

Attachment 2

Attachment 3

Attachment 4

Attachment 5

Attachment 6

Attachment 7

Attachment 8

Attachment 9

Attachment 10

Attachment 11

Attachment 12

Attachment 13

Attachment 14

Attachment 15

Attachment 16

Attachment 17

Attachment 18

Attachment 19

Attachment 20

Attachment 21

Attachment 22

Attachment 23

Attachment 24

Attachment 25

Attachment 26

Attachment 27

Attachment 28

Attachment 29

Attachment 30

Attachment 31

Attachment 32

Attachment 33