BACKGROUND
the study of human genome variation over the last several years has expanded from single nucleotide polymorphisms  to structural variants  of which copy number variants  constitute a distinct class  <cit> .
cnvs are unbalanced svs that are present in polymorphic copy number states, and include deletions, duplications, or combinations thereof.
recently, a comprehensive population-wide cnv map of the human genome was published  <cit> .
the authors reported on a large set of common cnvs  > 5%) greater than  <dig> kilobase  in length and genotyped  <dig>  of them on  <dig> individuals with ancestry in europe, africa or asia .
in addition, a recent report from the  <dig> genome project  <cit>  included about  <dig>  deletion polymorphisms  in the same hapmap populations from european, african and asian ancestry, and genotypes were inferred for most of the variants.
further mills et al.
characterized the majority of these deletions at sequence level resolution  <cit> , and proposed a detailed map of sv hotspots formed by common mechanisms.
cnvs have now been implicated with multiple common human diseases  including crohn’s disease   <cit> , osteoporosis   <cit> , body mass index  <cit> , and decreased susceptibility to hiv   <cit> .
cnvs have also been shown to be associated with high risk of autism  <cit> , schizophrenia  <cit>  and cancer  <cit> .
the effect of snps on gene expression was recently described through extensive expression quantitative trait locus  analysis in european and african populations  <cit> .
in  <dig>  stranger et al.
<cit>  performed transcripts association analysis both with snps and with cnvs, and highlighted similarities and differences for the two types of variation.
recently, two independent studies focused on the characterization of the potential impact of cnvs on gene transcripts by systematically querying cnv copy number states and paired gene expression levels obtained by rna sequencing on a common set of  <dig> individuals.
our group  <cit>  queried more than  <dig>  cnvs and showed that short cnvs  and gains are more likely to have functional impact with respect to larger cnvs and deletions, respectively.
schlattl et al.
<cit>  observed that cnvs exhibit a stronger correlation with expression than nearby snps and suggested a frequent causal role of cnvs in expression quantitative trait loci.
several human polymorphism studies, mainly based on snps, have been undertaken to identify human genome regions that are under selection .
the integrated map of positive selection from akey et al.
<cit>  suggested that positive selection targets encompass ~14% of the human genome and ~23% of all ucsc refseq genes.
in the context of snps, the signature of positive selection includes a high proportion of function-altering mutations, site frequency spectrum with high frequency of the derived allele and low genetic diversity , different allele frequency between populations, and long haplotypes  <cit> .
in the context of cnvs, positive selection could be suggested by allele frequencies that significantly differ between populations and by linkage disequilibrium with snps under positive selection  <cit> .
alternatively, population differentiation could be a result of genetic drift, however the latter is more significant for small populations  <cit> .
for instance, the salivary amylase gene  copy number varies considerably across populations and correlates with dietary starch prevalence, which supports the hypothesis of positive selection acting on amy <dig> copy number in high-starch diet populations  <cit> .
another example is the complex evolutionary history of the polymorphic ugt2b <dig> gene that shows high population differentiation  <cit>  .
recently, rhoxf <dig> has been reported as fast-evolving homeobox gene in primates with rapid evolution and copy number changes driven by darwinian positive selection acting on the male reproductive system  <cit> .
based on the observation of close genetic distances at the nucleotide level between human and chimpanzee, the hypothesis that regulatory mutations account for the majority of biological differences was proposed as early as  <dig>  <cit> .
more recent work suggests that gene expression levels can also serve as targets of selection  <cit> .
in the past, few studies focused on individual genes’ substitution rates  suggested that different types of selection act on genes depending on their position in the pathway  <cit> .
here, we posit more broadly that pathways can be subjects for selection.
existing models of pathway evolution, such as the horowitz retrograde model  <cit> , the chemistry-driven patchwork model  <cit> , and others  consider major modifications to a pathway chain, such as recruitment of a new enzyme  or a whole pathway duplication, that eventually lead to creation of a new pathway.
in addition to these major changes, we argue that more subtle ones might play an important role.
specifically, we hypothesize that fixed changes in gene product concentration levels - resulting from the changes in gene transcription levels - act as small adjustments while we assume that the number of pathway nodes and the pathway structure remain unchanged.
in this study, we performed an analysis of  <dig> well-characterized biological pathways seeking to determine how cnvs mapping to pathway genes impact transcript levels and how this effect differs across populations.
our results suggest that cnvs may modulate subtle changes in pathways at specific nodes, which may become fixed in certain populations.
we propose a model of pathway evolution where population differences are tuned at finite nodal points.
we discuss here the role of cnvs as potential modulators of biological pathways in human genome evolution.
methods
cnvs and hapmap genotype data
the complete set of  <dig>  cnv coordinates from  <cit>  was considered in the pathways enrichment/depletion analysis.
based upon the availability of high-resolution data  <cit>  at single sample level, we considered genotype calls for  <dig>  cnvs and then used this annotated set for population differentiation analysis.
the sample set included  <dig> ceu ,  <dig> yri ,  <dig> jpt , and  <dig> chb  from the international hapmap consortium .
throughout this study, japanese and chinese individuals are collectively referred to as asian .
all the data were downloaded from conrad et al.
<cit> .
classification of cnvs as gains and losses was taken exactly as inferred in  <cit> .
gene annotation and pathway information
refseq gene annotation information was downloaded from the university of california-santa cruz  web browser  <cit>  as ncbi build  <dig> .
the complete list of pathways and corresponding genes was compiled from kyoto encyclopedia of genes and genomes   <cit>  and biocarta .
f-statistics and cnv frequencies
population differentiation was evaluated using the f-statistics  <cit>  that ranges between  <dig>  and  <dig> .
f-statistics was calculated for each population pair by considering each cn genotype as an allele for diallelic cnvs as in mccarroll et al.
<cit> .
the f-statistics  was evaluated as follows: fst =  / ht ; ht =  <dig> - ∑ ti2; ti = +  ) / ; hs = · nx +  · ny ) / , where xi and yi are the population frequencies of allelic copy number i  in population x and y, respectively, nx and ny denote the number of individuals in population x and y, and ti is a weighted average of xi and yi.
this approach ignores the phase of the haplotype.
f-statistics cutoffs corresponding to the 5% and 1% tails of the distributions built for each population pair were considered.
specifically, the following values were calculated and applied: fst= <dig>   and fst= <dig>   for ceu-yri; fst= <dig>   and fst= <dig>  for  for yri-asn; fst= <dig>   and fst= <dig>   for ceu-asn.
beside the formal assessment of population differentiation by means of the f-statistics, we considered the frequency of polymorphisms, referred to as cnv frequency, as the sum of the frequencies of all cn states that differ from the major cn state.
cnv frequencies are utilized in the generation of frequency heatmaps.
for consistency, we selected the ceu major cn state for each cnv.
this choice does not affect the results as we focus on differences rather than absolute values.
size-dependent cnv enrichment/depletion analysis
given a set of cnvs and a set of genes  categories), the approach estimates if the number of observed cnv-gene overlaps is significantly higher  or lower  than what it would be expected by chance.to evaluate enrichment or depletion of pathways for the presence of overlapping cnvs, we developed a size-dependent enrichment method that takes into account the sizes of genes and cnvs and we empirically simulated the null distribution through permutations.
the approach overcomes the limitation of commonly applied statistical tests that do not consider either cnv or gene size information.
in reality, gene families widely differ in terms of gene sizes , respectively).
a schematic representation of this problem is presented in figure 1a where the blue  circles represent genes overlapping  with cnvs.
the following steps describe the method: first, for each chromosome, we compose a “pseudo transcriptome” made of concatenated genes  with flanking regions .
for simplicity, for genes with multiple isoforms coordinates corresponding to the widest size are considered.
we then define a ‘gene area’ as the area that spans the gene and upstream and downstream flanks.
second, we refine the cnv set by retaining only cnvs that overlap with gene areas.
in case of partial overlap, the cnv segment inside the gene area is retained .
third, the null distribution of enrichment is constructed by permuting this new set of truncated cnvs on the transcriptome by maintaining the number and size of cnvs per chromosome.
last, we define simple rules to count cnv-gene overlap as illustrated in figure 1c; the gene area is entirely within a cnv or the gene spans at least 50% of a cnv.
finally, for each pathway, the p-value is calculated by comparing the observed number of cnv-gene overlaps with the empirical null distribution constructed by the permutations.
it is important to mention few aspects that were not considered in the implementation of the null distribution.
first, the utilized cnvs genomic distribution information is limited by the use of array-based platforms.
second, paralogous genes, likely relevant in the context of segmental duplications where cnvs are significantly enriched, were not considered.
more refined approaches should eventually address these limitations.
the size dependent enrichment test executable and source code is downloadable at: http://demichelislab.unitn.it/tools/sizeenrichmenttest.tar.gz.
graphical representation of pathway related to cnv frequencies across populations
for each pathway, we considered the subset of genes overlapping with cnvs based on the gene areas  and represented it with the corresponding set of cnv-gene pairs.
based upon the systematic analyses of short- and long-range effect of cnvs on gene transcript levels  <cit> , we considered  <dig> kb and  <dig> mb flanks and composed two cnv-gene frequency maps.
the  <dig> kb flanks map will capture the majority of the significant associations between cnvs and genes  and the  <dig> mb flanks map will include longer distance effects.
we used  <dig> kb flanks for size-dependent enrichment/depletion analysis, and  <dig> mb flanks for cnv association analysis with gene expression levels.
for population differentiation analysis at the level of individual genes and of pathways we used  <dig> kb flanks to link genes and cnvs.
for population differentiation analysis of genes and pathways linked to “functional” cnvs we used  <dig> mb flanks.
one cnv can overlap multiple genes and vice-versa one gene can encompass  several cnvs.
therefore the number of cnv-gene pairs is not necessarily bound by the number of cnvs.
using  <dig> kb flanks we counted  <dig>  genes paired with cnvs and using  <dig> mb flanks we counted  <dig>  genes paired with cnvs.
heatmaps were used to graphically represent the cnv polymorphism frequency in each population, with cnv-gene pairs in rows, and populations in columns .
hierarchical clustering  was applied to identify groups of similar cnv-gene frequencies across the three populations.
gene expression data and association analysis between cnvs and gene expression levels
rna-seq data for ceu and  <dig> yri hapmap individuals were downloaded from  <cit>  and  <cit> .
sequencing data were originally generated using the illumina analyzer ii with 36-base and  <dig> or 46-base pairs, respectively.
yri individual data were downloaded from http://eqtl.uchicago.edu, ceu individual raw data were obtained from arrayexpress under accession numbers e-mtab- <dig> and preprocessed applying rseqtools  <cit> .
association analysis between gene expression data and cnvs was performed independently for ceu and yri individuals by cis analysis applying  <dig> mb flanks to each variant accordingly to  <cit> .
dosage effect  and allelic effect of transcript levels versus the copy number states were tested  <cit> .
multiple hypothesis testing correction was evaluated by calculating the false discovery rate  on the p-values; 10% and 5% thresholds were applied.
we will refer to cnv-gene pair with significant association between cnv states and gene expression levels as functional cnv-gene pair.
mechanisms of cnv formation
mechanisms of cnv formation were inferred as described in  <cit>  and are divided in four major classes: variable number of tandem repeats , non-allelic homologous recombination , transposable elements insertion  and non-homologous recombination .
repeatmasker annotation tracks were downloaded from ucsc browser to infer vntr and tei mechanisms.
two-sample test of proportions was applied to assess the significant differences in proportions of mechanism formation classes for cnvs differentiated in populations.
RESULTS
cnvs enrichment/depletion analysis at the level of pathways
we analyzed all annotated pathways from kegg  <cit>  and biocarta   for cnv enrichment.
for each pathway, we identified a set of cnv-gene pairs based on genomic co-location and applied size-dependent enrichment/depletion analysis .
the list of significantly enriched pathways is given in table  <dig>  and the distribution of enriched kegg functional classes and biocarta categories are presented in figure  <dig>  despite the general belief that cnvs are enriched mostly in extracellular and signaling pathways  <cit> , our results indicate significant enrichment in key metabolic pathway classes, such as carbohydrate , xenobiotic  and glycan  metabolisms .
among the metabolic pathways that are enriched in cnvs are the pathways from amino-acid  and energy  metabolisms.
* column fst lists the highest fst value for each population pair.
for the full list of pathways with cnv-gene pair annotation and corresponding fst values, see additional file 2: table s <dig>
enriched signaling pathways include keratinocyte differentiation  and phospholipase c d <dig> in phospholipid associated cell signaling  .
the list of depleted pathways and distribution of kegg and biocarta functional classes is given in additional file 2: table s <dig> and additional file 3: figure s <dig>  distribution of the depleted pathway classes is similar with that of all pathway classes in both databases.
focusing on the pathways enriched for the presence of cnvs, we compared the overall proportions of four cnv formation mechanisms and found evidence for enrichment of nhr  and depletion of tei  formation classes.
with respect to cnv types, gains or losses, depletion was observed for losses  and enrichment for gains .
population differentiation at the level of individual genes and of pathways
population differentiation analysis was performed for each cnv and each population pair-wise, by evaluating the f-statistics   <cit>  .
overall  <dig> differentiated cnv-gene pairs were detected at 5% significance level and  <dig> differentiated cnv-gene pairs were detected at 1% significance level, including previously reported genes  <cit> .
the cnv-gene pairs with the highest differentiation at 1% significance level are presented in table  <dig>  the distribution of cnv frequency differences for each pair of populations is presented in additional file 3: figure s <dig>  the data show that frequency differences are higher between yoruban and the two other populations than between europeans and asians, supporting the observation in  <cit> .
the same trend is observed in the distribution of allele frequency differences for snps  <cit> .
* the index refers to the original study where the gene has already been reported for population differentiation.
<dig> <cit> .
<dig> <cit> .
gene expression levels are also subject to selection in the human genome  <cit> .
at the level of individual genes, the selection acting on a cnv could be linked to transcription levels associated with the cnv states.
at the level of pathways, signatures of selection can be traced through the selection acting on individual genes constituting the pathway.
to investigate the extent of population differentiation on a pathway basis, we constructed pathway frequency maps  that combine polymorphism frequencies of each cnv-gene pair in each population and are visualized using a heatmap .
examples of pathway frequency maps are depicted in figure 3a and 3c.
focusing on the apoptosis pathway, we observed that akt <dig> and ikbkb related cnvs exhibit higher frequency of polymorphisms in the yoruban population compared to the european and asian populations.
tsta <dig>  b4galnt <dig> and sord cnv-gene pairs show significant population differentiation in the map for fructose and mannose metabolism, and ugt2b <dig>  ugt2b <dig>  ugt2a <dig> cnv-gene pairs are seen as highly differentiated in the map for androgen and estrogen metabolism c.
the analysis of all pathway cnv-gene frequency maps  suggested a characteristic pattern.
generally, a small proportion  of cnv-gene per pathway show strong population differentiation  while the remaining cnv-gene pairs are weakly differentiated or undifferentiated as we can observe for apoptosis, fructose and mannose metabolism and nf-kb signaling pathway .
a total of  <dig> pathways exhibited population differentiation in the sense that at least one cnv-gene pair is differentiated , among them  <dig> pathways showed the highest population differentiation .
this list includes the purine metabolism , the starch and sucrose metabolism , the pentose and glucuronate interconversions , and other highly differentiated metabolic pathways .
association of cnvs with gene expression at the level of genes and pathways
the proximity between a cnv and a gene can be considered an indicator for transcription regulation effect.
however, this only occurs for a fraction of variants.
to better assess the potential regulator effect of cnvs on gene transcripts, we considered rna sequencing gene expression data from  <dig> ceu  <cit>  and  <dig> yri  <cit>  individuals  as in  <cit> .
overall, we detected significant cnv gene expression association  involving  <dig> functional cnv-gene pairs in ceu of which  <dig> could be mapped to  <dig> pathways, and  <dig> functional cnv-gene pairs in yri of which  <dig> genes could be mapped to  <dig> pathways.
focusing on the set of functional cnv-gene pairs, we performed cnv enrichment/depletion analysis  for the set of  <dig> pathways in ceu and  <dig> pathways in yri  and found a significant enrichment of functional cnv-gene pairs for one pathway, the atrazine degradation pathway, in ceu and for  <dig>  mostly metabolic, pathways in yri.
the latter includes glutathione metabolism, starch and sucrose metabolism, cyanoamino acid metabolism, androgen and estrogen metabolism and others .
then we focused on functional cnv-gene pairs, which were characterized as significantly differentiated among populations.
out of  <dig>  functional cnv-gene pairs,  <dig>  showed significant population differentiation , and thus, are potential candidates for positive selection.
of the genes linked to functional cnvs that showed population differentiation ,  <dig> genes  are involved in  <dig> known pathways : the androgen and estrogen metabolism, the metabolism of xenobiotics by cytochrome p <dig>  the pentose and glucuronate interconversions, the porphyrin and chlorophyll metabolism, the starch and sucrose metabolism, the cell communication, the basal transcription factors, and the atrazine degradation.
the highest population differentiation  among functional cnvs was detected for  <dig> cnv-gene pairs: one is for the gene fam128a, and the other for the gene taf4b from the basal transcription factors pathway.
illustrations of significant associations between cnv and gene expression are presented for the gene taf4b belonging to the basal transcription factors pathway from yri data  and for the gene apobec3b belonging to the atrazine degradation pathway from ceu data .
taf4b encodes a subunit of transcription initiation factor that has been shown to function as co-activator of genes from nf-kb signaling pathway  <cit> .
functional cnv for taf4b is  <dig> kb large and located  <dig> kb upstream the gene of  <dig> kb.
observed cn states are  <dig>   <dig>  and  <dig>  and the expression is the highest for cn equal to  <dig>  cnv population frequencies are 89% in yri, 15% in ceu, and 42% in asn.
the gene apobec3b is a gene from cytidine deaminase family performing c to u rna-editing.
it has been known as an antiviral factor that can act against retroviruses, such as hiv  <cit> , and it has been also found to act as an inhibitor of l <dig> and alu retrotranspositon  <cit> .
functional cnvs for apobec3b is  <dig> kb large and completely removes the gene which is  <dig> kb in size.
population frequencies are 13% for ceu, 7% for yri, and 54% for asn.
illustrations of associations, though not statistically significant, between cnv and gene expression for genes from the mitogen-activated protein kinase  signaling pathway are presented in additional file 3: figure s <dig>  associations were detected for two genes located consecutively in the pathway chain: first gene is from the cacng family of gamma subunits important for regulating calcium channels  <cit> , and the second belongs to rasgrp family that activates map kinase cascade.
discussion
in this study, we performed a comprehensive size-dependent analysis of the impact of common cnvs on  <dig> biological pathways across human populations, with the idea that abundance of genetic variation in pathways can be indicative of ongoing evolution.
we first showed that cnv enriched pathways include not only signaling pathways and pathways involved in extracellular biological processes  <cit> , but also metabolic pathways such as amino-acid, carbohydrate, energy and glycan metabolisms.
the glycan biosynthesis pathway is annotated as originated in vertebrates in a recent phylogenetic study of metabolic pathways in the context of evolution  <cit>  and in the orangutan genome paper  <cit> .
some of the cnv-gene pairs we report are potential candidates for individual studies intended to investigate pathway-related dysfunctions and metabolic diseases  <cit>  .
in order to characterize the extent of human variation at the level of each pathway, we compiled a comprehensive list of population differentiated cnv-gene pairs.
as expected, the majority of pathways exhibiting population differentiation, belong to the signaling class, including calcium, mapk, toll-like receptor signaling, and others.
pathways from the metabolic class include purine metabolism, arginine and proline metabolism, and starch and sucrose metabolism with the well documented example of amy <dig> polymorphism and its relation to starch diet  <cit> .
these findings provide evidence that signaling pathways and pathways involved in extracellular activities are less conserved and are potentially under the influence of positive or adaptive selection  <cit> .
comprehensive characterization of the influence of negative selection is mandatory in future studies.
highly differentiated cnvs overlapping or in proximity of genes indicate recent evolutionary events, and emphasize the importance to improve our understanding of selection forces that shape the observed population differentiation  <cit> .
alternatively genetic drift events  <cit>  might be responsible for the observed population differentiation, however their effect is considered to be significant in small populations.
pathway based analysis revealed cnv-gene pairs with intermediate fst in addition to highly differentiated cnvs.
we reasoned that this might correspond to ‘cnvs genetic hitchhiking’ similar to what was suggested by barreiro  <cit>  in the context of snps.
recent studies also indicated the pervasiveness of negative, or purifying selection, acting on cnvs.
conrad et al.
<cit>  reported on purifying selection acting on exonic and intronic deletions.
other studies reported on variants under negative selection  <cit> .
the methods used for detection of purifying selection were based on the nucleotide level resolution, and restricted solely to the genes.
where our study focuses on the existence of significantly differentiated pattern, it is relevant to highlight that appropriate distinction between positive and purifying selection acting on cnvs is an important challenge that requires extensive future work.
we then investigated pathway enrichment considering the set of cnvs deemed functional through the association between the corresponding transcript levels and the copy number states and focused on those which also show differentiation across populations.
despite the fact that the transcript analysis is limited due to the sample size, the differences in the sequencing depth, and the sample type , we reasoned that our integrated analysis would benefit from this additional layer of information.
making a simplifying assumption, one can consider that selection acting on transcript levels ultimately influences the gene product leading to population specific selective advantage, as it was shown for ugt2b <dig> and amy <dig> in asians and europeans, respectively.
similarly, the variant associated with apobec3b  and present at different frequencies across populations may suggest a population specific antiviral effect, and in particular, population specific effect on alu-l <dig> retrotranspositional activity.
however, one has to consider that the detected potential “functional” cnvs are not necessarily direct causal variants and might simply be linked to the causal variants.
throughout the study some well-known cancer-related pathways were detected as enriched in population differentiated cnvs .
even though most of the detected genes have been extensively studied in relation to cancer, their differential effect across populations, including differential disease susceptibility, is still to be investigated.
the gene taf4b, whose expression levels are significantly associated with a population differentiated cnv  is a known co-activator of nf-kb genes further stimulating nf-kb transcriptional activity  <cit> .
it was shown that functional consequences of nf-kb signaling pathway is determined by nf-kb oscillation dynamics  <cit>  and the number, period and amplitude of nf-kb oscillations are regulated, via a negative feedback loop, by transcription levels of ikbα.
in a similar way, transcription levels of taf4b can regulate dynamics of nf-kb genes, providing different functional outcomes for the pathway.
mapk signaling pathway is another example of how the tuning of enzyme concentration can affect signaling pathway.
we found population differentiated functional cnvs for cancg and rasgrp <dig> located upstream of the mapk pathway chain .
the analysis of the model for this pathway identified existence of two different dynamical regimes, and depending on parameters, the system can switch from single-state bi-stability to oscillations  <cit> .
this means that the mapk signaling network can act as a switch and as a clock, and the altered  element concentration levels can initiate such transition or prevent cycling due to a shift in the threshold positions.
the effect of gene concentration level changes on the cellular phenotype and the concept of genetic balance has been addressed in  <cit> , together with the study of dosage-sensitive genes documented in genomes of various species, including yeast and human, and often encoding transcription regulators, signal transduction elements and binding factors  <cit> .
however, we recognized that the complexity of possible downstream changes relies on the non-linear dynamics of biochemical reactions possibly leading to non-proportional change in the concentration of the final component  <cit> .
here we provide insights into human pathways enriched for population differentiated functionally active cnvs.
under the assumption that a pathway chain remains intact as a whole , we hypothesize that evolutionary selected changes in transcription level of some genes constituting the pathway “tune” the pathway into a more favorable state for homeostasis, a process we refer to as the ‘tuning effect’.
last, we suggest that new pathways can stem as long-term potential outcome of the proposed tuning effect .
CONCLUSIONS
upon the characterization of functional cnvs and the concomitant population differentiation of the same variants suggestive of positive selection in different populations, it is challenging to discover the real effect of these changes on a pathway chain and to study the regulatory mechanisms in the cell that control the changes in gene concentration levels.
the picture becomes more complicated by acknowledging the multiplicity due to the existence of co-factors concurring to gene regulation, to the presence of other sources of variations, like epigenetic events, and by gene regulation compensatory effects.
our analysis may help to reveal pathway nodes, which have undergone changes  in gene concentrations, or, in other words, pathways that have been tuned.
further studies are required to understand the impact of these and other changes on pathway structure and human diversity.
competing interests
the authors declare that they have no competing interests.
authors’ contributions
mp carried out computational analyses, implemented the size-dependent enrichment analysis and drafted the manuscript.
sb carried out computational and statistical analysis and drafted the manuscript.
og participated in the interpretation of the results.
mar participated in the design of the study and interpretation of the results.
fd conceived the study, participated in its design, coordination and implementation and drafted the manuscript.
all authors read and approved the final manuscript.
supplementary material
additional file 1
supplementary material.
click here for file
additional file 2
table s <dig>  size-dependent enrichment analysis results for all the pathways considered in the study.
click here for file
additional file 3
figure s <dig>  distribution of depleted kegg and biocarta pathway classes.
a. distribution of depleted kegg classes.
original distribution of pathway classes in kegg database is given on the left.
distribution of the kegg classes depleted in cnvs that were obtained with the size-dependent enrichment analysis is given on the right.
b. distribution of the depleted biocarta categories.
original distribution of pathway categories in biocarta database is given on the left.
distribution of the biocarta categories depleted in cnvs that were obtained with the size-dependent enrichment analysis is given on the right.
figure s <dig>  histogram of cnv frequency differences in three population pairs: ceu-yri, ceu-asn and asn_yri.
cnv frequency is a frequency of polymorphism and calculated as described in the methods.
frequency differences are given in absolute values.
figure s <dig>  cnv-gene frequency heatmaps for  <dig> pathways.
heatmaps are constructed for cnv-gene pairs with  <dig> kb flanks.
figure s <dig>  rearrangements around sord gene area in human and chimpanzee.
figure shows mauve block alignment of four homologous regions in human and chimpanzee.
first two regions are extracted from human reference genome : chr15: <dig> ,000- <dig> , <dig> and chr15: <dig> ,000- <dig> , <dig>  and the second two regions are taken from chimpanzee genomes : chr15: <dig> ,000- <dig> , <dig> and chr15: <dig> ,000- <dig> , <dig>  the region of ~  <dig> kb that includes the gene sord underwent inverse duplication before the split of human and chimpanzee.
the active copy of the gene sord is encoded on the plus strand and is shown in the orange color.
the copy of the gene sord on the minus strand, that most likely became a pseudogene, is shown in the light orange.
cnv resulted from a loss of a region in human genome from the inverted copy of the gene sord .
analysis of the repeatmasker annotation revealed that in the chimpanzee, the l <dig> element  is located right next to the cnv boundary.
in the corresponding region in human, we see that the same l1pa <dig> element was truncated from 50% of length to 15%, and the alu element  was inserted just at the location of cnv.
however, only 62% of alu length remained in the sequence.
the transposable elements activity can also be seen in the promoter area of the gene sord.
the remnants of retrovirus  are present in the promoter region of three copies of sord except the active human copy .
also, full length l <dig> element  that was most likely inserted in the retrovirus, is observed upstream of the active copy of the gene sord, and the truncated copy of this element  remained upstream of the sord pseudogene copy.
figure s <dig>  mapk signaling pathway.
a. cnv-gene frequency heatmap for mapk signaling pathway.
rows correspond to cnv-gene pairs and columns correspond to three hapmap populations: yri, ceu and asn.
the values of the heatmap are cnv polymorphism frequency .
b. schematic representation of a fragment of mapk signaling pathway .
highlighted in orange are the gene families, cacn and rasgrp, whose genes have cnvs with significant gene expression associations  and evidence for population differentiation.
examples of gene variant associations for three genes from cacn family, cacng <dig>  cacng <dig> and cacng <dig>  are given in separate boxes.
example of gene variant association for rasgrp family is given for rasgrp <dig> gene.
figure s <dig>  tuning effect of pathway evolution.
different color and shape correspond to different enzymes.
increase in the concentration level of one enzyme  can induce changes in the concentration levels of the linked enzymes .
in the process of evolution, it can lead to the recruitment of enzymes that perform better functions, and as a result, create a new pathway.
click here for file
additional file 4
table s <dig>  complete list of  <dig> cnvs with fst calculated for each population pair.
click here for file
additional file 5
table s <dig>  list of pathways enriched for population differentiated cnv-gene pairs.
click here for file
additional file 6
table s <dig>  list of unique and shared pathways for venn diagram.
click here for file
additional file 7
table s <dig>  list of unique and shared genes for venn diagram of figure 4b.
click here for file
additional file 8
table s <dig>  enrichment analysis for the functional cnvs for ceu and yri.
click here for file
additional file 9
table s <dig>  list of cnvs that showed significant association with gene expression levels annotated with population differentiation statistics, pathways information and annotation on overlap with enhancers.
click here for file
acknowledgement
the authors would like to thank alex root, andrea sboner, and alessandro romanel for fruitful comments on this study and to acknowledge the department of defense new investigator award  and the starr cancer consortium for funding.

 the study of human genome variation over the last several years has expanded from single nucleotide polymorphisms to structural variants of which copy number variants constitute a distinct class. cnvs are unbalanced svs that are present in polymorphic copy number states, and include deletions, duplications, or combinations thereof. recently, a comprehensive population-wide cnv map of the human genome was published. the authors reported on a large set of common cnvs > 5%) greater than kilobase in length and genotyped of them on individuals with ancestry in europe, africa or asia. in addition, a recent report from the genome project included about deletion polymorphisms in the same hapmap populations from european, african and asian ancestry, and genotypes were inferred for most of the variants. further mills et al. characterized the majority of these deletions at sequence level resolution, and proposed a detailed map of sv hotspots formed by common mechanisms. cnvs have now been implicated with multiple common human diseases including crohn’s disease, osteoporosis, body mass index, and decreased susceptibility to hiv. cnvs have also been shown to be associated with high risk of autism, schizophrenia and cancer. the effect of snps on gene expression was recently described through extensive expression quantitative trait locus analysis in european and african populations. in stranger et al. performed transcripts association analysis both with snps and with cnvs, and highlighted similarities and differences for the two types of variation. recently, two independent studies focused on the characterization of the potential impact of cnvs on gene transcripts by systematically querying cnv copy number states and paired gene expression levels obtained by rna sequencing on a common set of individuals. our group queried more than cnvs and showed that short cnvs and gains are more likely to have functional impact with respect to larger cnvs and deletions, respectively. schlattl et al. observed that cnvs exhibit a stronger correlation with expression than nearby snps and suggested a frequent causal role of cnvs in expression quantitative trait loci. several human polymorphism studies, mainly based on snps, have been undertaken to identify human genome regions that are under selection. the integrated map of positive selection from akey et al. suggested that positive selection targets encompass ~14% of the human genome and ~23% of all ucsc refseq genes. in the context of snps, the signature of positive selection includes a high proportion of function-altering mutations, site frequency spectrum with high frequency of the derived allele and low genetic diversity, different allele frequency between populations, and long haplotypes. in the context of cnvs, positive selection could be suggested by allele frequencies that significantly differ between populations and by linkage disequilibrium with snps under positive selection. alternatively, population differentiation could be a result of genetic drift, however the latter is more significant for small populations. for instance, the salivary amylase gene copy number varies considerably across populations and correlates with dietary starch prevalence, which supports the hypothesis of positive selection acting on amy copy number in high-starch diet populations. another example is the complex evolutionary history of the polymorphic ugt2b gene that shows high population differentiation. recently, rhoxf has been reported as fast-evolving homeobox gene in primates with rapid evolution and copy number changes driven by darwinian positive selection acting on the male reproductive system. based on the observation of close genetic distances at the nucleotide level between human and chimpanzee, the hypothesis that regulatory mutations account for the majority of biological differences was proposed as early as. more recent work suggests that gene expression levels can also serve as targets of selection. in the past, few studies focused on individual genes’ substitution rates suggested that different types of selection act on genes depending on their position in the pathway. here, we posit more broadly that pathways can be subjects for selection. existing models of pathway evolution, such as the horowitz retrograde model, the chemistry-driven patchwork model, and others consider major modifications to a pathway chain, such as recruitment of a new enzyme or a whole pathway duplication, that eventually lead to creation of a new pathway. in addition to these major changes, we argue that more subtle ones might play an important role. specifically, we hypothesize that fixed changes in gene product concentration levels - resulting from the changes in gene transcription levels - act as small adjustments while we assume that the number of pathway nodes and the pathway structure remain unchanged. in this study, we performed an analysis of well-characterized biological pathways seeking to determine how cnvs mapping to pathway genes impact transcript levels and how this effect differs across populations. our results suggest that cnvs may modulate subtle changes in pathways at specific nodes, which may become fixed in certain populations. we propose a model of pathway evolution where population differences are tuned at finite nodal points. we discuss here the role of cnvs as potential modulators of biological pathways in human genome evolution. methods cnvs and hapmap genotype data the complete set of cnv coordinates from was considered in the pathways enrichment/depletion analysis. based upon the availability of high-resolution data at single sample level, we considered genotype calls for cnvs and then used this annotated set for population differentiation analysis. the sample set included ceu, yri, jpt, and chb from the international hapmap consortium. throughout this study, japanese and chinese individuals are collectively referred to as asian. all the data were downloaded from conrad et al.. classification of cnvs as gains and losses was taken exactly as inferred in. gene annotation and pathway information refseq gene annotation information was downloaded from the university of california-santa cruz web browser as ncbi build. the complete list of pathways and corresponding genes was compiled from kyoto encyclopedia of genes and genomes and biocarta. f-statistics and cnv frequencies population differentiation was evaluated using the f-statistics that ranges between and. f-statistics was calculated for each population pair by considering each cn genotype as an allele for diallelic cnvs as in mccarroll et al.. the f-statistics was evaluated as follows: fst = / ht; ht = - ∑ ti2; ti = + ) /; hs = · nx + · ny ) /, where xi and yi are the population frequencies of allelic copy number i in population x and y, respectively, nx and ny denote the number of individuals in population x and y, and ti is a weighted average of xi and yi. this approach ignores the phase of the haplotype. f-statistics cutoffs corresponding to the 5% and 1% tails of the distributions built for each population pair were considered. specifically, the following values were calculated and applied: fst= and fst= for ceu-yri; fst= and fst= for for yri-asn; fst= and fst= for ceu-asn. beside the formal assessment of population differentiation by means of the f-statistics, we considered the frequency of polymorphisms, referred to as cnv frequency, as the sum of the frequencies of all cn states that differ from the major cn state. cnv frequencies are utilized in the generation of frequency heatmaps. for consistency, we selected the ceu major cn state for each cnv. this choice does not affect the results as we focus on differences rather than absolute values. size-dependent cnv enrichment/depletion analysis given a set of cnvs and a set of genes categories), the approach estimates if the number of observed cnv-gene overlaps is significantly higher or lower than what it would be expected by chance.to evaluate enrichment or depletion of pathways for the presence of overlapping cnvs, we developed a size-dependent enrichment method that takes into account the sizes of genes and cnvs and we empirically simulated the null distribution through permutations. the approach overcomes the limitation of commonly applied statistical tests that do not consider either cnv or gene size information. in reality, gene families widely differ in terms of gene sizes, respectively). a schematic representation of this problem is presented in figure 1a where the blue circles represent genes overlapping with cnvs. the following steps describe the method: first, for each chromosome, we compose a “pseudo transcriptome” made of concatenated genes with flanking regions. for simplicity, for genes with multiple isoforms coordinates corresponding to the widest size are considered. we then define a ‘gene area’ as the area that spans the gene and upstream and downstream flanks. second, we refine the cnv set by retaining only cnvs that overlap with gene areas. in case of partial overlap, the cnv segment inside the gene area is retained. third, the null distribution of enrichment is constructed by permuting this new set of truncated cnvs on the transcriptome by maintaining the number and size of cnvs per chromosome. last, we define simple rules to count cnv-gene overlap as illustrated in figure 1c; the gene area is entirely within a cnv or the gene spans at least 50% of a cnv. finally, for each pathway, the p-value is calculated by comparing the observed number of cnv-gene overlaps with the empirical null distribution constructed by the permutations. it is important to mention few aspects that were not considered in the implementation of the null distribution. first, the utilized cnvs genomic distribution information is limited by the use of array-based platforms. second, paralogous genes, likely relevant in the context of segmental duplications where cnvs are significantly enriched, were not considered. more refined approaches should eventually address these limitations. the size dependent enrichment test executable and source code is downloadable at: http://demichelislab.unitn.it/tools/sizeenrichmenttest.tar.gz. graphical representation of pathway related to cnv frequencies across populations for each pathway, we considered the subset of genes overlapping with cnvs based on the gene areas and represented it with the corresponding set of cnv-gene pairs. based upon the systematic analyses of short- and long-range effect of cnvs on gene transcript levels, we considered kb and mb flanks and composed two cnv-gene frequency maps. the kb flanks map will capture the majority of the significant associations between cnvs and genes and the mb flanks map will include longer distance effects. we used kb flanks for size-dependent enrichment/depletion analysis, and mb flanks for cnv association analysis with gene expression levels. for population differentiation analysis at the level of individual genes and of pathways we used kb flanks to link genes and cnvs. for population differentiation analysis of genes and pathways linked to “functional” cnvs we used mb flanks. one cnv can overlap multiple genes and vice-versa one gene can encompass several cnvs. therefore the number of cnv-gene pairs is not necessarily bound by the number of cnvs. using kb flanks we counted genes paired with cnvs and using mb flanks we counted genes paired with cnvs. heatmaps were used to graphically represent the cnv polymorphism frequency in each population, with cnv-gene pairs in rows, and populations in columns. hierarchical clustering was applied to identify groups of similar cnv-gene frequencies across the three populations. gene expression data and association analysis between cnvs and gene expression levels rna-seq data for ceu and yri hapmap individuals were downloaded from and. sequencing data were originally generated using the illumina analyzer ii with 36-base and or 46-base pairs, respectively. yri individual data were downloaded from http://eqtl.uchicago.edu, ceu individual raw data were obtained from arrayexpress under accession numbers e-mtab- and preprocessed applying rseqtools. association analysis between gene expression data and cnvs was performed independently for ceu and yri individuals by cis analysis applying mb flanks to each variant accordingly to. dosage effect and allelic effect of transcript levels versus the copy number states were tested. multiple hypothesis testing correction was evaluated by calculating the false discovery rate on the p-values; 10% and 5% thresholds were applied. we will refer to cnv-gene pair with significant association between cnv states and gene expression levels as functional cnv-gene pair. mechanisms of cnv formation mechanisms of cnv formation were inferred as described in and are divided in four major classes: variable number of tandem repeats, non-allelic homologous recombination, transposable elements insertion and non-homologous recombination. repeatmasker annotation tracks were downloaded from ucsc browser to infer vntr and tei mechanisms. two-sample test of proportions was applied to assess the significant differences in proportions of mechanism formation classes for cnvs differentiated in populations.  cnvs enrichment/depletion analysis at the level of pathways we analyzed all annotated pathways from kegg and biocarta for cnv enrichment. for each pathway, we identified a set of cnv-gene pairs based on genomic co-location and applied size-dependent enrichment/depletion analysis. the list of significantly enriched pathways is given in table and the distribution of enriched kegg functional classes and biocarta categories are presented in figure despite the general belief that cnvs are enriched mostly in extracellular and signaling pathways, our results indicate significant enrichment in key metabolic pathway classes, such as carbohydrate, xenobiotic and glycan metabolisms. among the metabolic pathways that are enriched in cnvs are the pathways from amino-acid and energy metabolisms. * column fst lists the highest fst value for each population pair. for the full list of pathways with cnv-gene pair annotation and corresponding fst values, see additional file 2: table s enriched signaling pathways include keratinocyte differentiation and phospholipase c d in phospholipid associated cell signaling. the list of depleted pathways and distribution of kegg and biocarta functional classes is given in additional file 2: table s and additional file 3: figure s distribution of the depleted pathway classes is similar with that of all pathway classes in both databases. focusing on the pathways enriched for the presence of cnvs, we compared the overall proportions of four cnv formation mechanisms and found evidence for enrichment of nhr and depletion of tei formation classes. with respect to cnv types, gains or losses, depletion was observed for losses and enrichment for gains. population differentiation at the level of individual genes and of pathways population differentiation analysis was performed for each cnv and each population pair-wise, by evaluating the f-statistics. overall differentiated cnv-gene pairs were detected at 5% significance level and differentiated cnv-gene pairs were detected at 1% significance level, including previously reported genes. the cnv-gene pairs with the highest differentiation at 1% significance level are presented in table the distribution of cnv frequency differences for each pair of populations is presented in additional file 3: figure s the data show that frequency differences are higher between yoruban and the two other populations than between europeans and asians, supporting the observation in. the same trend is observed in the distribution of allele frequency differences for snps. * the index refers to the original study where the gene has already been reported for population differentiation... gene expression levels are also subject to selection in the human genome. at the level of individual genes, the selection acting on a cnv could be linked to transcription levels associated with the cnv states. at the level of pathways, signatures of selection can be traced through the selection acting on individual genes constituting the pathway. to investigate the extent of population differentiation on a pathway basis, we constructed pathway frequency maps that combine polymorphism frequencies of each cnv-gene pair in each population and are visualized using a heatmap. examples of pathway frequency maps are depicted in figure 3a and 3c. focusing on the apoptosis pathway, we observed that akt and ikbkb related cnvs exhibit higher frequency of polymorphisms in the yoruban population compared to the european and asian populations. tsta b4galnt and sord cnv-gene pairs show significant population differentiation in the map for fructose and mannose metabolism, and ugt2b ugt2b ugt2a cnv-gene pairs are seen as highly differentiated in the map for androgen and estrogen metabolism c. the analysis of all pathway cnv-gene frequency maps suggested a characteristic pattern. generally, a small proportion of cnv-gene per pathway show strong population differentiation while the remaining cnv-gene pairs are weakly differentiated or undifferentiated as we can observe for apoptosis, fructose and mannose metabolism and nf-kb signaling pathway. a total of pathways exhibited population differentiation in the sense that at least one cnv-gene pair is differentiated, among them pathways showed the highest population differentiation. this list includes the purine metabolism, the starch and sucrose metabolism, the pentose and glucuronate interconversions, and other highly differentiated metabolic pathways. association of cnvs with gene expression at the level of genes and pathways the proximity between a cnv and a gene can be considered an indicator for transcription regulation effect. however, this only occurs for a fraction of variants. to better assess the potential regulator effect of cnvs on gene transcripts, we considered rna sequencing gene expression data from ceu and yri individuals as in. overall, we detected significant cnv gene expression association involving functional cnv-gene pairs in ceu of which could be mapped to pathways, and functional cnv-gene pairs in yri of which genes could be mapped to pathways. focusing on the set of functional cnv-gene pairs, we performed cnv enrichment/depletion analysis for the set of pathways in ceu and pathways in yri and found a significant enrichment of functional cnv-gene pairs for one pathway, the atrazine degradation pathway, in ceu and for mostly metabolic, pathways in yri. the latter includes glutathione metabolism, starch and sucrose metabolism, cyanoamino acid metabolism, androgen and estrogen metabolism and others. then we focused on functional cnv-gene pairs, which were characterized as significantly differentiated among populations. out of functional cnv-gene pairs, showed significant population differentiation, and thus, are potential candidates for positive selection. of the genes linked to functional cnvs that showed population differentiation, genes are involved in known pathways: the androgen and estrogen metabolism, the metabolism of xenobiotics by cytochrome p the pentose and glucuronate interconversions, the porphyrin and chlorophyll metabolism, the starch and sucrose metabolism, the cell communication, the basal transcription factors, and the atrazine degradation. the highest population differentiation among functional cnvs was detected for cnv-gene pairs: one is for the gene fam128a, and the other for the gene taf4b from the basal transcription factors pathway. illustrations of significant associations between cnv and gene expression are presented for the gene taf4b belonging to the basal transcription factors pathway from yri data and for the gene apobec3b belonging to the atrazine degradation pathway from ceu data. taf4b encodes a subunit of transcription initiation factor that has been shown to function as co-activator of genes from nf-kb signaling pathway. functional cnv for taf4b is kb large and located kb upstream the gene of kb. observed cn states are and and the expression is the highest for cn equal to cnv population frequencies are 89% in yri, 15% in ceu, and 42% in asn. the gene apobec3b is a gene from cytidine deaminase family performing c to u rna-editing. it has been known as an antiviral factor that can act against retroviruses, such as hiv, and it has been also found to act as an inhibitor of l and alu retrotranspositon. functional cnvs for apobec3b is kb large and completely removes the gene which is kb in size. population frequencies are 13% for ceu, 7% for yri, and 54% for asn. illustrations of associations, though not statistically significant, between cnv and gene expression for genes from the mitogen-activated protein kinase signaling pathway are presented in additional file 3: figure s associations were detected for two genes located consecutively in the pathway chain: first gene is from the cacng family of gamma subunits important for regulating calcium channels, and the second belongs to rasgrp family that activates map kinase cascade. discussion in this study, we performed a comprehensive size-dependent analysis of the impact of common cnvs on biological pathways across human populations, with the idea that abundance of genetic variation in pathways can be indicative of ongoing evolution. we first showed that cnv enriched pathways include not only signaling pathways and pathways involved in extracellular biological processes, but also metabolic pathways such as amino-acid, carbohydrate, energy and glycan metabolisms. the glycan biosynthesis pathway is annotated as originated in vertebrates in a recent phylogenetic study of metabolic pathways in the context of evolution and in the orangutan genome paper. some of the cnv-gene pairs we report are potential candidates for individual studies intended to investigate pathway-related dysfunctions and metabolic diseases. in order to characterize the extent of human variation at the level of each pathway, we compiled a comprehensive list of population differentiated cnv-gene pairs. as expected, the majority of pathways exhibiting population differentiation, belong to the signaling class, including calcium, mapk, toll-like receptor signaling, and others. pathways from the metabolic class include purine metabolism, arginine and proline metabolism, and starch and sucrose metabolism with the well documented example of amy polymorphism and its relation to starch diet. these findings provide evidence that signaling pathways and pathways involved in extracellular activities are less conserved and are potentially under the influence of positive or adaptive selection. comprehensive characterization of the influence of negative selection is mandatory in future studies. highly differentiated cnvs overlapping or in proximity of genes indicate recent evolutionary events, and emphasize the importance to improve our understanding of selection forces that shape the observed population differentiation. alternatively genetic drift events might be responsible for the observed population differentiation, however their effect is considered to be significant in small populations. pathway based analysis revealed cnv-gene pairs with intermediate fst in addition to highly differentiated cnvs. we reasoned that this might correspond to ‘cnvs genetic hitchhiking’ similar to what was suggested by barreiro in the context of snps. recent studies also indicated the pervasiveness of negative, or purifying selection, acting on cnvs. conrad et al. reported on purifying selection acting on exonic and intronic deletions. other studies reported on variants under negative selection. the methods used for detection of purifying selection were based on the nucleotide level resolution, and restricted solely to the genes. where our study focuses on the existence of significantly differentiated pattern, it is relevant to highlight that appropriate distinction between positive and purifying selection acting on cnvs is an important challenge that requires extensive future work. we then investigated pathway enrichment considering the set of cnvs deemed functional through the association between the corresponding transcript levels and the copy number states and focused on those which also show differentiation across populations. despite the fact that the transcript analysis is limited due to the sample size, the differences in the sequencing depth, and the sample type, we reasoned that our integrated analysis would benefit from this additional layer of information. making a simplifying assumption, one can consider that selection acting on transcript levels ultimately influences the gene product leading to population specific selective advantage, as it was shown for ugt2b and amy in asians and europeans, respectively. similarly, the variant associated with apobec3b and present at different frequencies across populations may suggest a population specific antiviral effect, and in particular, population specific effect on alu-l retrotranspositional activity. however, one has to consider that the detected potential “functional” cnvs are not necessarily direct causal variants and might simply be linked to the causal variants. throughout the study some well-known cancer-related pathways were detected as enriched in population differentiated cnvs. even though most of the detected genes have been extensively studied in relation to cancer, their differential effect across populations, including differential disease susceptibility, is still to be investigated. the gene taf4b, whose expression levels are significantly associated with a population differentiated cnv is a known co-activator of nf-kb genes further stimulating nf-kb transcriptional activity. it was shown that functional consequences of nf-kb signaling pathway is determined by nf-kb oscillation dynamics and the number, period and amplitude of nf-kb oscillations are regulated, via a negative feedback loop, by transcription levels of ikbα. in a similar way, transcription levels of taf4b can regulate dynamics of nf-kb genes, providing different functional outcomes for the pathway. mapk signaling pathway is another example of how the tuning of enzyme concentration can affect signaling pathway. we found population differentiated functional cnvs for cancg and rasgrp located upstream of the mapk pathway chain. the analysis of the model for this pathway identified existence of two different dynamical regimes, and depending on parameters, the system can switch from single-state bi-stability to oscillations. this means that the mapk signaling network can act as a switch and as a clock, and the altered element concentration levels can initiate such transition or prevent cycling due to a shift in the threshold positions. the effect of gene concentration level changes on the cellular phenotype and the concept of genetic balance has been addressed in, together with the study of dosage-sensitive genes documented in genomes of various species, including yeast and human, and often encoding transcription regulators, signal transduction elements and binding factors. however, we recognized that the complexity of possible downstream changes relies on the non-linear dynamics of biochemical reactions possibly leading to non-proportional change in the concentration of the final component. here we provide insights into human pathways enriched for population differentiated functionally active cnvs. under the assumption that a pathway chain remains intact as a whole, we hypothesize that evolutionary selected changes in transcription level of some genes constituting the pathway “tune” the pathway into a more favorable state for homeostasis, a process we refer to as the ‘tuning effect’. last, we suggest that new pathways can stem as long-term potential outcome of the proposed tuning effect.  upon the characterization of functional cnvs and the concomitant population differentiation of the same variants suggestive of positive selection in different populations, it is challenging to discover the real effect of these changes on a pathway chain and to study the regulatory mechanisms in the cell that control the changes in gene concentration levels. the picture becomes more complicated by acknowledging the multiplicity due to the existence of co-factors concurring to gene regulation, to the presence of other sources of variations, like epigenetic events, and by gene regulation compensatory effects. our analysis may help to reveal pathway nodes, which have undergone changes in gene concentrations, or, in other words, pathways that have been tuned. further studies are required to understand the impact of these and other changes on pathway structure and human diversity. competing interests the authors declare that they have no competing interests. authors’ contributions mp carried out computational analyses, implemented the size-dependent enrichment analysis and drafted the manuscript. sb carried out computational and statistical analysis and drafted the manuscript. og participated in the interpretation of the results. mar participated in the design of the study and interpretation of the results. fd conceived the study, participated in its design, coordination and implementation and drafted the manuscript. all authors read and approved the final manuscript. supplementary material additional file 1 supplementary material. click here for file additional file 2 table s size-dependent enrichment analysis results for all the pathways considered in the study. click here for file additional file 3 figure s distribution of depleted kegg and biocarta pathway classes. a. distribution of depleted kegg classes. original distribution of pathway classes in kegg database is given on the left. distribution of the kegg classes depleted in cnvs that were obtained with the size-dependent enrichment analysis is given on the right. b. distribution of the depleted biocarta categories. original distribution of pathway categories in biocarta database is given on the left. distribution of the biocarta categories depleted in cnvs that were obtained with the size-dependent enrichment analysis is given on the right. figure s histogram of cnv frequency differences in three population pairs: ceu-yri, ceu-asn and asn_yri. cnv frequency is a frequency of polymorphism and calculated as described in the methods. frequency differences are given in absolute values. figure s cnv-gene frequency heatmaps for pathways. heatmaps are constructed for cnv-gene pairs with kb flanks. figure s rearrangements around sord gene area in human and chimpanzee. figure shows mauve block alignment of four homologous regions in human and chimpanzee. first two regions are extracted from human reference genome: chr15:,000-, and chr15:,000-, and the second two regions are taken from chimpanzee genomes: chr15:,000-, and chr15:,000-, the region of ~ kb that includes the gene sord underwent inverse duplication before the split of human and chimpanzee. the active copy of the gene sord is encoded on the plus strand and is shown in the orange color. the copy of the gene sord on the minus strand, that most likely became a pseudogene, is shown in the light orange. cnv resulted from a loss of a region in human genome from the inverted copy of the gene sord. analysis of the repeatmasker annotation revealed that in the chimpanzee, the l element is located right next to the cnv boundary. in the corresponding region in human, we see that the same l1pa element was truncated from 50% of length to 15%, and the alu element was inserted just at the location of cnv. however, only 62% of alu length remained in the sequence. the transposable elements activity can also be seen in the promoter area of the gene sord. the remnants of retrovirus are present in the promoter region of three copies of sord except the active human copy. also, full length l element that was most likely inserted in the retrovirus, is observed upstream of the active copy of the gene sord, and the truncated copy of this element remained upstream of the sord pseudogene copy. figure s mapk signaling pathway. a. cnv-gene frequency heatmap for mapk signaling pathway. rows correspond to cnv-gene pairs and columns correspond to three hapmap populations: yri, ceu and asn. the values of the heatmap are cnv polymorphism frequency. b. schematic representation of a fragment of mapk signaling pathway. highlighted in orange are the gene families, cacn and rasgrp, whose genes have cnvs with significant gene expression associations and evidence for population differentiation. examples of gene variant associations for three genes from cacn family, cacng cacng and cacng are given in separate boxes. example of gene variant association for rasgrp family is given for rasgrp gene. figure s tuning effect of pathway evolution. different color and shape correspond to different enzymes. increase in the concentration level of one enzyme can induce changes in the concentration levels of the linked enzymes. in the process of evolution, it can lead to the recruitment of enzymes that perform better functions, and as a result, create a new pathway. click here for file additional file 4 table s complete list of cnvs with fst calculated for each population pair. click here for file additional file 5 table s list of pathways enriched for population differentiated cnv-gene pairs. click here for file additional file 6 table s list of unique and shared pathways for venn diagram. click here for file additional file 7 table s list of unique and shared genes for venn diagram of figure 4b. click here for file additional file 8 table s enrichment analysis for the functional cnvs for ceu and yri. click here for file additional file 9 table s list of cnvs that showed significant association with gene expression levels annotated with population differentiation statistics, pathways information and annotation on overlap with enhancers. click here for file acknowledgement the authors would like to thank alex root, andrea sboner, and alessandro romanel for fruitful comments on this study and to acknowledge the department of defense new investigator award and the starr cancer consortium for funding.
we implemented an enrichment analysis of pathways accounting for cnvs and genes sizes and detected significant enrichment not only in signal transduction and extracellular biological processes, but also in metabolism pathways. the majority of these cnv-gene pairs belong to signal transduction pathways and 6% of the cnv-gene pairs show statistical association between the copy number states and the transcript levels. the analysis suggested possible examples of positive selection within individual populations including nf-kb, mapk signaling pathways, and alu/l retrotransposition factors. altogether, our results suggest that constitutional cnvs may modulate subtle pathway changes through specific pathway enzymes, which may become fixed in some populations. inherited copy number variants can modulate the expression levels of individual genes. however, little is known about how cnvs alter biological pathways and how this varies across different populations. upon the estimation of cnv population differentiation, we evaluated that 22% of the pathways contain at least one gene that is proximal to a cnv that shows significant population differentiation. to trace potential evolutionary changes of well-described biological pathways, we jointly queried the genomes and the transcriptomes of a collection of individuals with caucasian, asian or yoruban descent combining high-resolution array and sequencing data.
 inherited copy number variants can modulate the expression levels of individual genes. however, little is known about how cnvs alter biological pathways and how this varies across different populations. to trace potential evolutionary changes of well-described biological pathways, we jointly queried the genomes and the transcriptomes of a collection of individuals with caucasian, asian or yoruban descent combining high-resolution array and sequencing data.  we implemented an enrichment analysis of pathways accounting for cnvs and genes sizes and detected significant enrichment not only in signal transduction and extracellular biological processes, but also in metabolism pathways. upon the estimation of cnv population differentiation, we evaluated that 22% of the pathways contain at least one gene that is proximal to a cnv that shows significant population differentiation. the majority of these cnv-gene pairs belong to signal transduction pathways and 6% of the cnv-gene pairs show statistical association between the copy number states and the transcript levels.  the analysis suggested possible examples of positive selection within individual populations including nf-kb, mapk signaling pathways, and alu/l retrotransposition factors. altogether, our results suggest that constitutional cnvs may modulate subtle pathway changes through specific pathway enzymes, which may become fixed in some populations. cnvspathwayspathway evolutionpopulation geneticseqtl 