Retrogene survival is not impacted by linkage relationships

Data collection

Retrogene data was collected from http://yeti.amu.edu.pl/retrogenedb/ (accessed 3/21/21) (Kabza, Ciomborowska & Makałowska, 2014). Regulatory network data was collected from RegNetwork for humans, and JASPAR and REDfly for D. melanogaster (accessed 1/17/21) (Liu et al., 2015; Mathelier et al., 2013; Rivera et al., 2018). Assembly size, chromosome sizes, and genomic coordinates for each of the genes were collected from NCBI (human: GRCh38.p13, D. melanogaster: FB release 6). After removing retrocopies with missing information (IDs, protein-coding status, etc), a total of 4,426 retrocopies were found for humans, with 106 retrogenes. Retrocopies were paired with 1,431 parental genes with 809 unique network partners. In D. melanogaster, 82 retrocopies were found, 81 of which were retrogenes. These retrocopies were subsequently coupled with 64 parental genes with 109 unique network partners. Genetic distances were estimated using recombination maps for intrachromosomal parent-retrogene pairs; for interchromosomal pairs, distance was set to a default value of 0.5. (Kong et al., 2002; Rezvoy et al., 2007).

There are a number of significant differences between the human and D. melanogaster dataset. As mentioned, the human data contains many more pseudoretrogenes than the D. melanogaster data; this difference is likely an artifact of data collection and prior research directions rather than an indication of the true rate of pseudogenization. Most of the data on human pseudogenes originates from Ciomborowska et al. (2013), which extensively analyzes retrocopy content in humans; similar studies have not been conducted in D. melanogaster. Additionally, the density of regulatory networks differs between humans and D. melanogaster, with the latter exhibiting much denser and well-connected networks. Whether this is biological truth is a different question; however, for this study, we assume that the network data we have is an unbiased sample from the true relationships.

Hypothesis testing

All analysis was done using R 3.9.1 in the RStudio 1.2.1335 IDE (R Core Team, 2017; RStudio Team, 2020). Packages used include dplyr, rentrez, chromPlot, GenomicFeatures, MareyMap and reshape2 (Lawrence et al., 2013; Rezvoy et al., 2007; Verdugo & Orostica, 2019; Wickham, 2007; Wickham, 2016; Wickham et al., 2021). The reference distributions were generated using a Monte Carlo permutation method, with the null assumption that retrogenes are i.i.d. from a uniform random distribution over the genome. We had originally considered constructing the reference distributions by sampling retropseudogenes from our data; however, we decided against this for two reasons: firstly, the distribution of retropseudogenes observed in our data may not accurately reflect the true distribution of retrocopies at the time of retroduplication, and secondly, the distribution of retropseudogenes appears to be approximately uniform across the genome regardless (Fig. S1).

Samples of random retrogenes are simulated by drawing genomic coordinates from the above distribution. Relevant statistics (described later in the Results) are then calculated. For each of our tests, this process was iterated 1,000 times to produce an approximate reference distribution for the relevant statistic. Test statistics were then calculated using the empirical data and compared against the reference distribution using a two-tailed test (α = 0.05). All code and data used in this study are available via the author’s GitHub repository (https://github.com/johnathanlo/retrogenes).

Testing data quality

Before proceeding with analysis, the quality of the dataset was assessed by searching for evidence of the well-documented out-of-the-X effect. As previously mentioned, this is a pattern in which more retrogenes originate from the X chromosome than expected through random chance. The statistic of interest is

$$X_1={\displaystyle \{\mathrm{r}\mathrm{e}\mathrm{t}\mathrm{r}\mathrm{o}\mathrm{g}\mathrm{e}\mathrm{n}\mathrm{e}\mathrm{s}\mathrm{f}\mathrm{r}\mathrm{o}\mathrm{m}\mathrm{t}\mathrm{h}\mathrm{e}\mathrm{X}\mathrm{c}\mathrm{h}\mathrm{r}\mathrm{o}\mathrm{m}\mathrm{o}\mathrm{s}\mathrm{o}\mathrm{m}\mathrm{e}\}\cap \{\mathrm{r}\mathrm{e}\mathrm{t}\mathrm{r}\mathrm{o}\mathrm{g}\mathrm{e}\mathrm{n}\mathrm{e}\mathrm{s}\mathrm{n}\mathrm{o}\mathrm{t}\mathrm{o}\mathrm{n}\mathrm{t}\mathrm{h}\mathrm{e}\mathrm{X}\mathrm{c}\mathrm{h}\mathrm{r}\mathrm{o}\mathrm{m}\mathrm{o}\mathrm{s}\mathrm{o}\mathrm{m}\mathrm{e}\}.}$$

Reference distribution was constructed by simulating random pairs of retrogenes and parental genes, then counting the number of pairs with a parental gene located on the X and a retrogene located on a different chromosome. In each realization of the Monte Carlo simulation, 106 retrogene/parental gene pairs were generated for humans, and 81 pairs for D. melanogaster, corresponding to each of the retrogenes in their respective datasets. This procedure was repeated 1,000 times to generate the null distribution. After calculating the test statistic from our data, significant results were found in both species, with p ∼ 0 (Fig. S1). A key assumption behind our testing procedure is the assumption that retrocopies and parental genes are distributed randomly in the genome.

Retrogene distribution with respect to parental genes in humans

We wanted to test the linkage relationships for 106 retrogenes in humans and 81 retrogenes in D. melanogaster. To determine whether or not retrogenes are distributed independently of parental location, we analyzed the average genetic distance between a retrogene and its parental gene. The out-of-the-X effect was controlled for by removing all such parental genes/retrogenes from both datasets. Since our model posits a uniform random distribution of retrogenes independent of parental gene location, this should not bias our results on a subset of strictly autosomal parental genes.

The statistic is

$$X_2={\displaystyle \frac{\sum _{allb}dist(a,b)}{\frac{}{}}n}$$where a is the parental gene, b is the retrogene, and n is the sample size (83 and 55 for humans and D. melanogaster respectively). The distance function is computed by applying recombination maps to the genomic coordinates of the parent/retrogene pair, with a default value of 0.5 for pairs on different chromosomes. The reference distribution was constructed as before by simulating a set of random parental genes and a corresponding set of random retrogenes and calculating the above statistic for 10,000 realizations. The test statistics were calculated from the data. No significant results were observed for either human or D. melanogaster (Fig. 2).

The above statistic weights each retrogene equally; however, there is heterogeneity in the number of retrogenes generated by each parental gene. To account for any bias resulting from this, we additionally test whether or not the average retrogene for each parental gene is more or less genetically distant than expected; in other words, we weight distances by parental genes instead. For each parental gene, we may calculate the statistic

$$X_3={\displaystyle \frac{\sum dist(a,x)}{|C|}}$$where a is the parental gene and x ∈ C, where C is the set of retrocopies associated with a. X₃ is therefore a function of the random variable dist (a, x). We define the mean over the population of parental genes as

$$\mu ={\displaystyle E(X_3).}$$

To obtain a distribution over the sample mean, we simulate the statistic X₃ using Monte Carlo methods as before, then compute the estimate

$$\hat{\mu }={\displaystyle \frac{1}{n}\sum _{i=0}^n{x}_i}$$in the usual way, where x_i assumed to be i.i.d. instances of X₃ as defined above, and test at a significance level of α = 0.05. The results show that lack of significance persists regardless of whether genetic distance is weighted by retrogenes or by parental genes (Fig. 3). In other words, the average retrogene survival event is not influenced by linkage, and the average parental gene does not produce a distribution of retrogenes significantly different from random (null) expectation.

Retrogene distribution relative to network partners

A final test was performed to assess whether retrogenes are randomly distributed with respect to their nearest network partners. Since we hypothesized that the linkage relationships of the parent may affect the fixation patterns of retrogenes, it is reasonable to extend this influence to network partners of the parent. Retrogene fixation may be influenced by gene regulatory network topology, since close proximity of a gene to network partners can ensure that certain combinations of alleles are inherited together as members of a co-adapted gene complex. Hence, we examine the relationship of retrogenes to the nearest network partner of parental genes, including the parental gene itself. We define a statistic describing the minimum distance between retrogenes and network partners as

$$X_4={\displaystyle min\{dist(a,x):x\in C\}}$$where a is a retrogene, x is a network member/parental gene, C is the set of network partners of the parental gene and the parental gene, and dist(a, x) is the genetic distance between a and x as defined in the previous section; note dist (a, x) is a random variable, so X₄ is random. Then we define the mean over the population

$$\mu ={\displaystyle E(X_4)}$$and construct a distribution using Monte Carlo permutations as before. From the sample, we can then compute an estimate of the mean

$$\hat{\mu }={\displaystyle \frac{1}{n}\sum _{i=0}^n{x}_i}$$where we assume x_i i.i.d. samples from X₄ and test at a significance level of α = 0.05. This test finds no significant deviation of the sample statistic from the expectation under the null for either humans or D. melanogaster (Fig. 4). In other words, the distribution of retrogenes in the genome does not seem to be influenced by proximity to parental network partners.

Linkage patterns do not play significant role in retrogene fixation

None of our tests exhibit any significant association between changes in linkage and retrogene fixation in either humans or D. melanogaster. In fact, the results markedly conform to our null expectation of uniformly random distribution of retrogenes across the genome. In other words, given that a retrogene fixes, it does not have any tendency to be closer or further to either its parental gene or the network partners of its parental genes.

Linkage with nearest network partners does not affect retrogene fixation

The relationship between retrogenes and network partners was also investigated. We investigated the linkage relationship between retrogenes and their nearest network partner specifically. No significant correlations were observed in either humans or D. melanogaster. These findings do not exclude the possibility of a broader effect involving network partners. The nearest network partner is insufficient to fully characterize the topology of an entire regulatory network, and it remains plausible that the overall topology of the network is correlated with changes in linkage from parental gene to retrogene.

Concluding Remarks

This study provides theoretical background and a preliminary investigation into a novel hypothesis regarding retrogene evolution. While these results do not indicate a significant role for linkage in determining retrogene fixation, several limitations and confounding factors provide basis for further investigation. One confounding factor is that variation in selection on linkage (e.g., selection for increased linkage in some lineages vs selection for decreased linkage in others) may mask signals from detection by our methods. Indeed, we know that in mammals and Drosophila, the out-of-the-X pattern is a definitive example of selection against linkage, while other studies have demonstrated that proximity may sometimes be selected for to derive the benefits of nearby regulatory regions or open chromatin formations (Bai, Casola & Betrán, 2008; Loppin et al., 2005). Additionally, even though we find no consistent pattern, selection on linkage may still have been key to the survival of a subset of retrogenes. Our study assumes that retrogenes are i.i.d. with respect to genetic distances, a simplification that may not hold up in reality. To uncover effects of linkage on the level of individual genes would require in-depth functional and comparative analysis of suspected retrogene/parental gene pairs similar to work on the out-of-the-X retrogenes. The analysis is of course limited also by constraints related to available data; data from diverse taxa may be necessary to sufficiently illuminate the role of linkage in retrogene survival.

Finally, it remains interesting to ask whether or not selection on linkage can be mediated through retrogene fixation, much like the argument for the fixation of modifier alleles. The primary obstacle to such a study is finding a sample of genes that have experienced selection on their linkage relationships at some point in the past. With such a sample, we can ask whether or not they produce retrogenes at a greater than normal rate. Constructing such a sample with any certitude may seem like a daunting task, but one that may be amenable through experimental evolution techniques.

There are plausible reasons for why selection on linkage may play a role in retrogene fixation only rarely. Previous work has strongly emphasized neofunctionalization as an outcome of retrogene fixation, which would make such retrogenes less likely to interact with the parental gene or its network partners (Casola & Betrán, 2017; He & Zhang, 2005). Our work suggests exactly this: the observed lack of correlation between parental gene networks and retrogenes may indicate that retrogenes typically occupy different regulatory networks and fulfill different functions when compared to their parents. This phenomenon may also be diminished by limitations on the expression of new retrocopies. New retrocopies require expression to be selected for, and since they do not typically carry regulatory elements with them, they may not achieve consistent or appreciable levels of expression, which prevents selection from acting. The primary venue for new retrocopies to achieve high levels of expression is during promiscuous expression, as in the thymus and testes, or during haploid expression during spermatogenesis (Casola & Betrán, 2017; Raices, Otto & Vibranovski, 2019). Then, retrocopies may only be selected for if their expression conveys fitness benefits under these particular contexts, and retrocopies that do not provide any effect or benefit during spermatogenesis would thus be less likely to experience selection.

As genomic sequencing and analysis tools improve, analysis of retrogene evolutionary dynamics in other species will become feasible. Further investigation of the hypotheses presented here might be best served by data from deer and muntjac (Mudd et al., 2020; Yang et al., 1997). Their shared evolutionary history makes comparisons between findings between species feasible; additionally, these lineages have undergone multiple chromosome fusions and fissions in recent evolutionary time, which provide a backdrop of changes in linkage that make some form of selection on linkage nearly inevitable. Another potentially fruitful line of inquiry would be to ask if specific classes of parental genes and retrogenes are more influenced by selection on linkage than others. For example, recent work has found autosomal pairs of parental genes and retrogenes that exhibit complementarity of expression in the testes, which may appear to be a more likely scenario for selection on linkage (Casola & Betrán, 2017). Both of the above present reasonable directions for future work in this area.