Manual Gene Retention and Fractionation Bias Data for 37 Panels of Orthologous Sorghum-Maize1-Maize2 Chromosomal Segments

To study the post-tetraploidy fractionation process in detail, both a sequenced genome that is undergoing fractionation and an outgroup with a sequenced genome that diverged before the tetraploidy event is required. For the most recent tetraploidy in maize, which happened from 5 to 12 MYA [15],[16], sorghum is such an outgroup. Sorghum diverged from the maize lineage just before the tetraploidy event (Figure 1) [17]. Sorghum has been sequenced [16], and the first assembly of a maize genome has recently been published [13]. In addition to its phylogenetic position, the maize lineage tetraploidy possesses another characteristic that recommends it for the study of fractionation: it is relatively recent. The alpha tetraploidy in Arabidopsis, at 23–50 MYA, is older than the most recent (alpha) tetraploidy in maize. Even so, the maize alpha tetraploidy is known to be highly fractionated [10],[18]. Our primary research aim was to detail what happened to those orthologous (syntenic) genes shared by sorghum and maize following the maize tetraploidy.

Using the procedures described in Methods, we identified 37 orthologous regions between sorghum and the corresponding maize homeologs retained after the maize alpha tetraploidy event. From these regions, we found that of the 2,943 sorghum-maize (Sb-Zm) syntenic shared genes that we studied, 43% of them were retained as homeologous pairs in maize. Note that we count as present any significant fragment of gene. If the maize tetraploidy behaved as other known tetraploidies in plants and microbes, retained genes should be enriched in those encoding transcription factors, as reviewed [19]. Indeed, the frequency of genes encoding transcription factors was 4.3 times greater among the retained genes as compared to the fractionated genes. Figure 3A is a cartoon of a GEvo output screenshot of a 13-gene segment of one of the 37 orthologous regions (region Sb2) between sorghum and its two maize homeologs, as described in Methods. The GEvo comparative sequence alignment tool output generated the original blastn output detailed in Methods. Supplemental Information 1 (Dataset S1) gives our primary data as inferred from analyses like that shown in Figure 3A. Dataset S2 shows how any one sorghum chromosomal region is orthologous to two maize regions, generating information essential to construct our sorghum-maize1-maize2 regions.

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Figure 3. Determining which maize homeologous genes were deleted following tetraploidy.

(A) Cartoon of a GEvo blastn comparison graphic depicting a 13-gene stretch of sorghum with the two orthologous regions of maize. Sorghum nucleotides were masked except those in official gene models or genes shared between sorghum and maize (SI1); maize DNA was masked for repeat sequences. Blastn high scoring pairs (hits) are colored orange if they are Sb-Zm1 and purple if they are Sb-Zm2. Colored lines indicate orthology. The code B021BB2DDDB10B abbreviates these data, where B = both genes remain, 0 = gene missing in syntenous Zm positions (these are discarded), 1 = the gene on Zm1 remains alone because its homeolog was deleted, 2 = the gene on Zm2 remains alone because its homeolog was deleted, and D = local duplicates of an arbitrary mother gene leftmost in the cluster (D's are discarded). Therefore, the essential code for this 13-gene Sb stretch reduces to B21BB2B1B. The circle indicates a Zm1 gene that has some completely and some partially deleted exons; we noted these partially fractionated genes for further research, but we counted them as present (B). The brackets enclose clusters of tandemly duplicated genes in both Sb and Zm2. The arrow indicates a single gene in Sb hitting a reverse tandem duplication in maize; maize genes like this one were counted as present. Please use http://genomevolution.org/r/37e to reproduce on-the-fly the Sb2 region blast experiment, the region containing Sb02g030760-Sb02g030950 drawn above. (B) Using the same color code of the panel above, these are two exemplary small regions from the total of 37 regions comprising our syntenic dataset. The regions exemplifying fractionation bias (AL3–7) and no bias (AL3–8) are color-coded in such a way that the number of gene deletions suffered by Zm homeolog 1 versus Zm homeolog 2 is easy to count. A 1∶1 (p = 0.5) ratio of these deletions computes to not biased. p≤0.1 is weak support and ≤0.05 is strong support that these deletions significantly deviate from this 1∶1 null hypothesis.

https://doi.org/10.1371/journal.pbio.1000409.g003

One way to measure fractionation bias is to first assume that gene loss involves one gene independently from any other gene, and then count the number of gene losses (deletions) on one homeolog as compared to the number of gene losses on the other homeolog. If fractionation were unbiased, this ratio is expected to be 1∶1. Other measures of fractionation include total number of genes or base pairs in an orthologous stretch, but counting deletions of shared genes is most direct, so we present this first. Figure 3B shows two representative diagrams of our data for shorter regions (Dataset S3 contains all 37 such diagrams with bias statistics) and indicates that fractionation has been significantly biased in 68% of our regions. Using data from nine representative longer sorghum regions (Table 1), we conclude that the over-fractionated chromosomes have 2.3 times as many deleted genes as do the under-fractionated chromosomes.

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Table 1. The organization of fractionation and similar data based on the designation Zm-under-fractionated (u) homeolog and Zm over-fractionated homeolog (o).

https://doi.org/10.1371/journal.pbio.1000409.t001

We next asked, what was the average extent of gene loss? Most importantly, are deletion events longer than one gene? Figure 4A, B, C, and D monitors runs, or the sequential series of deleted genes, and Figure 4E monitors runs of retained genes. The experimental runs data plotted in Figure 4A, B, C, and E were compared to the 95% confidence interval around the median of Monte Carlo simulated data (Methods) based on the assumption that one gene is deleted at a time, and the chromosome was ligated before the next deletion, as such a mechanism would be predicted to work in nature. The distributions of Figures 4A, B, and C are all very similar: they differ only as to whether or not the over-fractionated or under-fractionated chromosome is evaluated or as to whether or not gene losses of 10 genes or greater were included. The most frequent run length in distributions Figure 4A, B, or C is one gene, followed by two genes, and so forth. If we recalculate expectations for distribution of Figure 4C using an evolutionary method that permits varying percentages of deletions of genes, the best fit is one-gene deletions 80%, two-gene deletions 15%, and three-gene deletions 5% (Figure 4D).

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Figure 4. Distribution of runs: Runs of lost genes, or runs of genes retained.

x-axis: length of the run, in genes. y-axis: number of runs. Black bars represent observed data. The white circles within the gray bars represent the median number of runs of that length observed from 1,000 simulations. The gray bars mark the limits within which the values of 95% from the simulations fell. Panels A–D refer to runs of deletions. Panels E and F refer to runs of genes retained from the maize lineage tetraploidy. (A) Observed and distributions of deletion runs in the over-fractionated homeologous regions including all deletion data, and simulated distributions assuming genes are lost solely through a 1 gene deletion mechanism. (B) Observed distributions of deletion runs in the over-fractionated homeologous regions with deletions longer than nine genes removed as likely segmental transpositions, and simulated distributions assuming genes are lost solely through a 1 gene deletion mechanism. (C) Observed and simulated distributions of deletions in the under-fractionated homeologous regions, assuming genes are lost solely through a 1 gene deletion mechanism. (D) Observed and simulated distributions of deletion runs in the over-fractionated homeologous regions, using a model where 80% of deletions are single gene, 15% remove two adjacent genes, and 5% of deletions remove three genes in a row, best fit ratio determined by a genetic algorithm. (E) Observed runs of genes retained in both homeologs and simulated distributions assuming genes are lost solely through a 1 gene deletion mechanism. (F) Observed runs of genes retained in both homeologs compared to the simulated expectation assuming 80% single-gene deletions, 15% two-gene, and 5% three-gene.

https://doi.org/10.1371/journal.pbio.1000409.g004

The possibility existed that longer deletion runs were not authentic deletions but were segmental translocations. Deletion runs consisting of 12 genes or more were found somewhere else in the genome. Deletions of between 11 and 6 genes were found elsewhere about 10% of the time, but identification was made more difficult because fractionation is expected to remove genes from any position in the genome, and sometimes is expected to leave behind fewer than the three syntenic genes needed for a positive identification. That is why deletions of 10 or more genes were removed from all distributions of Figure 4 except Figure 4A. There is a possibility, a possibility we evaluate, that the smaller deletions are also undetectable segmental translocations and not authentic deletions, but removal is essentially one gene at a time.

Clear Differentiation Between Deletion of a Segment of DNA and Translocation of That Segment to a New Position

We next asked if it were possible that, rather than being deleted, single genes observed as lost between orthologous maize and sorghum regions were instead transposed or translocated elsewhere in the genome, as we had observed for longer runs of genes. Large-scale single gene transposition has been documented in the eudicot order Brassicales [9], and cases have also been reported in the maize lineage [11]. To address the possibility that the majority of the fractionating gene loss we were observing was actually a result of whole-gene transposition, we attempted to identify potentially orthologous maize genes in a position-independent manner. Sorghum genes with known orthologs in rice were blasted against sorghum, rice, and maize genomes, as described in Methods. From the resulting data we found that genes identified as retained had a mean of 1.57 copies in the genome, with a median of 2 genes. It is expected that this number would be less than 2, as the manual annotators considered a gene to be retained if a significant fragment of it was still present, which included genes in the process of being removed by small deletions (as will be discussed). Genes identified as being fractionated in Dataset S1 were present at a mean of 1.17 copies in the genome, with a median of 1 gene. A few of these extra copies are likely maize-specific duplications, but others no doubt represent apparently deleted homeologs that have transposed to other locations in the genome. Nevertheless, these data provide strong evidence that while some apparently fractionated genes may have been lost via translocation (transposition to a new site), translocation is not the prevailing mechanism explaining our fractionation data. This conclusion does not imply that fragments of genes are not transposed around the genome, as is known to occur frequently via transposon-mediated gene capture [20]. Indeed, when we re-calibrated our search to find shorter stretches of high-identity sequence, we found many pieces of genes present at higher copy numbers elsewhere in the genome. Examination of a sample of these hits identified gene fragments, but no intact genes were found.

Genes Retained as Pairs from the Maize Lineage Tetraploidy Are Rarely Clustered Beyond Expectations

As shown in Figure 2, fractionation itself clusters retained genes. Figure 4E identifies runs of retained genes (Bs) and distributes them by run length and compares this to expectations based on deletions one gene at a time. Is this distribution more highly clustered than expected from fractionation alone? The mode is clearly one retained pair, as expected. Expectation intervals were generated assuming that deletions occurred one gene at a time. Although clustering of retained genes is not dramatic, runs of retained genes greater than 9 gene pairs are not expected at all; in total, there are 62 genes (out of 1,203, or 5%) in such longer runs ranging from 9 through 12 genes in length. When expectations are changed to be 80% single-gene deletions, 15% two-gene deletions, and 5% three-gene deletions (Figure 4F), the actual and expected are similar. Now there are only four unexpected runs greater than 9 genes in length. With the exception of these few longer runs, genes are retained approximately as expected based on 80%/15%/5% 1/2/3 gene deletion predictions.

Use of Whole-Gene Count Data to Evaluate Fractionation Bias Throughout the Genome

Table 1 focuses on nine longer representative homeologous regions of maize representing different sorghum chromosomes. The under-fractionated and over-fractionated homeologs in maize are identified in this table (Column 2). This over/under designation derives from the deletion bias data quantified in Table 1 and evaluated for significance in Column C. Table 2 shows these data for each of the nine representative regions individually (Column H). In this case, the numbers are less than 1 because the ratio is under-fractionated/over-fractionated, and the under-fractionated homeolog has fewer deletions. The homeolog with the fewest deletions contains the most genes, so another measure of fractionation is the number of genes on the under-fractionated/over-fractionated, where bias will now be indicated by ratios greater than 1. The fractionation bias ratios, using total gene data, for each of the nine representative regions are listed in Column L of Table 1.

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Table 2. Repeat sequence signature analyses of 16 Zm deletions within exons of eight genes.

https://doi.org/10.1371/journal.pbio.1000409.t002

To extrapolate bias in our manually annotated regions to more of the genome, we used the slope of syntenic lines in Zm-Zm dot plots (Dataset S4). A slope of 1 implies unbiased fractionation. A significant difference in the number of genes or base pairs between the two homeologous maize chromosomes alters that slope from 1, and this is what we observed. If the unit of Zm-under/Zm-over measurement is total number of genes annotated by maizegnome.org, the average slope value corresponds to a mean fractionation bias value of 1.5 (Table 1, Column M). If the units are in total base-pairs, the fractionation bias is 2.3 (Table 1, Column N). Again, the under/over direction of fractionation in both cases remains greater than 1, as expected, but the dot-plot analysis made it possible to examine considerably longer regions of paired homeologs, each anchored on the indicated sorghum region. Most important here is that the three measures of bias based on gene number (Columns L and M) or base pair length (N) are concordant with expectations based on the rigorous deletion bias data generated manually for our representative regions. Based on the concordance between bias in orthologous gene loss and base pair length, and given that 85% of the maize genome is composed of transposable elements [13], we conclude that homeologous regions that preferentially lose genes also lose intergenic, primarily transposon, DNA more frequently.

The Under- and Over-Fractionated Maize Homeologs Are Equally Diverged from Their Sorghum Ortholog at the Level of Nearly Neutral Base Pair Substitution (Ks)

Table 1, Column O, reports our measured ratio of Ks [Zm-under/Sb] to Ks [Zm-over/Sb] for a total of 1,772 Sb-Zm-Zm gene units in the five sorghum homeologous regions for which highly significant under/over-fractionation expectations existed (Table 1). We removed 16% of pairs with the most extreme Ks ratios, many of which represent misalignments or alignments to pseudogenes. Using the remaining data, we found no difference between the Ks values between sorghum and either of the two maize homeologous regions. We conclude that mutation by base substitution and mutation by short deletion are mechanically distinct and are targeted differently.

There Is No Obvious Correlation between an Over-Fractionated Chromosomal Arm and the Number of Map Units (% Recombination) in That Arm

Three of our representative regions are within pairs of homeologous chromosomal arms: Sb1 = Zm1S/9L (sorghum chromosome 1 = maize chromosome 1S and maize chromosome 9L), Sb3 = Zm3L/8L, and Sb6 = Zm2S/10L. The under-fractionated (longer) homeolog is the numerator. These are the only arm-arm exact homeologies in the maize genome; examination of syntenic Sb-Zm dot plots (like that in Dataset S2) made clear that segments of these arms are not present syntenically on any other chromosomes. The total map unit's length of these maize chromosomal arms is known, making it possible to directly compare the degree of fractionation within any given arm to the overall recombination frequency within that arm. Mapping data for maize inbred T232×x inbred CM37 generated the following data [21]: the proportion of map units for under-fractionated arm/over-fractionated arms are Zm1S/9L = 0.9, Zm3L/8L = 1.1, and Zm2S/10L = 1.9 (Table 1, Column N). Note that although Zm2S/10L has the largest difference in recombination frequency, it has the lowest fractionation bias of these three paired arms (Table 1, Column C). We conclude that there is no obvious correlation between biased fractionation and overall frequencies of reciprocal recombination during meiosis.

Quantitative Estimate of Methylation Status for Over- and Under-Fractionated Maize Homeologous Chromosomal Segments

Even before BAC sequencing was complete, one group [22] identified methylation domains of maize chromosome in shoot or root nuclei using McrBC restriction endonuclease, a treatment that degrades DNA between methylated half sites of the form m5C-N40–500-m5C. McrBC is non-specific for different types of methylation patterns. Using this crude measure of methylated regions (BAC start-stop) in maize shoot nuclei, we overlifted (translated the start-stop nucleotide designations) the data from BACs to pseudomolecules and found no correlation at all between the over-fractionated and under-fractionated homeologs (Table 1, Column M). Two representative regions were concordant, two were not concordant, and one region was vastly over-methylated on the under-fractionated homeolog.

Single Deletions within Exons Occur Primarily on the Over-Fractionated Chromosomes and Appear to Be Due to Illegitimate Recombination

Our methods for deciding whether or not a maize gene was retained (“B”) did not require that the entire coding sequence be present, but only a significant fragment. Because of this, our calculation of the number of whole sorghum-maize genes retained post-tetraploidy, about 40%, is surely an overestimate. If we were to assume that the process of fractionation is ongoing, we reasoned that some of our retained genes might have internal deletions whose flanking sequences might give us a clue as to the mechanism behind gene fractionation. By visual examination, we identified cases where a maize gene seemed to have a gap within an exon. To verify each fully flanked deletion, we extracted the sorghum exon sequence and used it as query for a blastn to rice, a grass that diverged from sorghum about 50 MYA [23], to sorghum itself, and to the two homeologous maize regions. We then studied each Os-Sb-Zm1-Zm2 blastn result using GEvo, our synteny visualization platform (in CoGe, Methods). We verified that eight genes, containing a total of 16 deletions, were fully flanked by conserved, known sequence. Table 2 gives the data for these fully flanked deletions. Figure 5A and B shows an exemplary GEvo graphic and the pertinent orthologous sequences of rice, sorghum, and the two maize homeologs. In two cases (Sb01g039030 and Sb09g023840, Table 2), the apparent gap was actually several short gaps within the homeologous flanking sequence. The gap size within these 16 deleted regions ranged from 5 bp to 178 bp, with a mean gap size of 25.9 bp. Bias for gaps is consistent with the fractionation bias found locally: in other words, when a gap is present, it is in the maize homeologous gene located on over-fractionated chromosome 93% (15/16) of the time (Table 2).

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Figure 5. Deletions and intra-chromosomal recombination in fractionated maize homeologs.

(A) A cartoon of a blastn output between orthologous genes in rice, sorghum, and two maize homeologs. The colored rectangles (blue, orange, purple) represent high-scoring sequence pairs (HSPs), or regions with high sequence similarity to each other. Sorghum is the reference sequence (Sb09). The top two panels show the syntenous exon sequence (highlighted in blue) between rice (rice chromosome 5, or Os05) and sorghum (sorghum chromosome 9, or Sb09); the second and fourth panels demonstrate the location of the deletion in maize homeolog Zm chromosome 8 (Zm08) to its sorghum ortholog (purple). As can be seen, the deletion is evident when the Zm08 sequence (circled) is compared to the orthologous sequence of Zm06 (orange) and rice (blue) when all HSPs are aligned on sorghum. The GEvo alignment output for these data may be found at http://genomevolution.org/r/3em. (B) A ClustalW alignment of the rice, sorghum, and maize homeologs from Figure 5A. The purple sequence in the unfractionated maize homeolog (Zm06) indicates the location of the direct repeat sequence that originally flanked the deletion in the fractionated homeolog (Zm08). (C) A diagram representing the mechanism of intra-chromosomal recombination, based on the flanking sequence highlighted in (B). Direct repeats come together to form a circle, which is then recombined away, leaving a solo repeat in its place.

https://doi.org/10.1371/journal.pbio.1000409.g005

As mentioned in the Introduction, deletions due to illegitimate recombination are often flanked by a short stretch of sequence that, before the deletion, had been a direct repeat [11]. In theory, such repeats facilitate ectopic, intra-chromosomal, reciprocal recombination (as drawn in Figure 5C) generating a circle and a solo copy of the original repeat sequence in place of the sequence deleted (the circle). Using ClustalW, we found such direct repeats flanking 10 of the 16 gaps in our study (Table 2). These repeats were between 3 and 24 bp in length; an example is given in Figure 5B. Notice how the repeats surrounding the gap in the fractionated homeolog are truncated in comparison to the repeat sequence within the whole homeolog: this is a typical footprint of intra-chromosomal recombination [11].

In an attempt to generalize our results from monocots (e.g. grasses) to eudicots (e.g. legumes), we found several such small deletions where the inferred precursor sequence was flanked by direct repeats, within retained duplicate genes of Glycine max (Gm: soybean, unpublished data) from the more recent of the two easily observable tetraploid events in the sorghum genome. Soybean has had two recent genome duplication events, the most recent one (alpha) having taken place between 14 and 3 MYA [14]. The close relative, Medicago trunculata, was used as the outgroup in order to infer the precursor gene sequence before deletion. We conclude that small deletions are involved in the fractionation of genes following ancient (successful) plant tetraploidies.

Plant Genomic Sequences

The sorghum sequence was Sbi1 assembly and Sbi1.4 annotation (Paterson et al. 2009 [23]) downloaded from Phytosome V4.0 http://www.phytozome.net/sorghum, last modified 3-25-08. The B73 maize genome sequence was obtained in the form of pseudomolecules in 3-09 (ftp://ftp.genome.arizona.edu/pub/fpc/maize/) and stored in our CoGe platform as database 8082: http://synteny.cnr.berkeley.edu/CoGe/OrganismView.pl?dsgid=8082; and with draft models annotations in 10-09 from maizegenome.org (http://ftp.maizesequence.org/release-4a.53/). The sequence of these two releases is identical. The draft annotated maize sequence will be called “4a.53” in the few instances where we use the official CDS models. The TIGR 5 Nipponbare rice assembly and annotation was downloaded onto our CoGe platform (http://synteny.cnr.berkeley.edu/CoGe/) before the MSU6 update file://localhost/ (http://rice.plantbiology.msu.edu/data_download.shtml), and was used in 2008 to generate the sorghum gene list used here; the differences between rice TIGR5 and MSU6 annotations are of no consequence to this project. The soybean and Medicago trunculata genomes were downloaded from Phytosome 4.0 http://www.phytozome.net/soybean.php and http://www.phytozome.net/medicago.php, respectively, in early 2009.

The Sorghum Gene List (Dataset S1)

Sorghum genes were Sbi1.4 to which we added many genes on the basis of orthology to rice Nipponbare, TIGR 5; the added genes included many with corresponding RNAs since these are absent in Sbi1.4. Dataset S1 uses the format Sbxgxxxxxx for Sbi1.4 genes and sorghum_chrmosomex_startx_stopx for genes we added based on Sb-Os orthology. The detailed syntenic alignment of the entire genomes of sorghum and rice was automated and frozen in September 2009 as the rice-sorghum CNS discovery Pipeline 1.0. Dataset S6 diagrams this pipeline and details each step. What is most important is that any sorghum gene we use in this analysis is shared syntenically between sorghum and least one of the two possible homeologous maize positions. Some of our added sorghum genes are shared with maize as orthologs; those that are not shared with maize were not studied.

After adding 10,585 new putative genes to the 34,003 official JGI (Joint Genome Institute) sorghum genes, the augmented sorghum genome was masked for any sequence repeated over 50 times, and then everything but exons or RNA-encoding sequence was additionally masked. This heavily masked sorghum genome was then used to query the maize genome.

We found a total of 37 orthologous regions between sorghum and the corresponding maize homeologs retained after the maize alpha tetraploidy event (Sb-Zm1-Zm2) in two ways, and both ways used applications available online in the CoGe comparative genomics platform (http://synteny.cnr.berkeley.edu/CoGe/). Central to our success was our ability to clearly visualize the locations of the many translocations and inversions that happened in both the sorghum and maize lineages. Knowing all breakpoints makes it clear that any single sorghum chromosome is orthologous to exactly two maize chromosome regions, even though many smaller segments are often involved (Dataset S2). To this end all of the 37 regions begin and end with at least one gene retained by both maize homeologs. In this way (Methods), a total of 4,461 sorghum genes (10% of the sorghum genome) were set up for manual evaluation. In order to define those genes that had an ortholog in maize, we condensed all members of locally duplicated arrays into one gene and discarded those 492 duplicates (11%), leaving one parent gene for each array. We also invalidated 74 genes that had annotation incongruencies, and then disregarded another 953 genes for which we failed to find maize orthologs. Each sorghum gene is given an evaluation code of “1” (has an ortholog in the first Zm homeologous region), “2” (has an ortholog in the second Zm homeologous region), or “B” (has an ortholog in both Zm homeologous regions). The designations “O,” “N,” and “D”: D = local duplicate; N = invalid data; O = no ortholog in Zm. In each of the 37 regions within Dataset S1, every sorghum gene has been annotated with one of these six symbols. A link (tinyurl.com, a URL abbreviation service, or genomevolution.org) is provided for each Sb-Zm1-Zm2 panel to facilitate the repetition of our research in the GEvo alignment graphic tool we used for research (Dataset S1). We were left with 2,943 orthologous Sb-Zm1-Zm2 genes spread over 37 Sb-Zm1-Zm2 regions.

Using Applications in the CoGe Toolbox, Including the Ks Values Provided in SynMap

CoGe is an integrated collection of maintained databases, algorithms and applications useful to compare complete genomes on demand [28],[29] and without which it would be difficult to perform our analyses in a reasonable amount of time. SynMap, within CoGe, is a dot plot application that implements the DAGchainer algorithm [30] to identify syntenic lines in two-dimensional arrays of blastn hits between two identical (to find homeologies) or different (to find orthologies) genomes. Each “dot” is a gene pair. The color of this dot can be portrayed to reflect Ks (synonymous base-pair substitution frequency), so syntenic lines of different ages have different colors (see Dataset S2). Clicking on any dot in SynMap anchors the GEvo sequence comparison tool and automatically generates a blastn alignment output. Each output (like a BLAST or LAGAN output) includes a graph, a link, and can be repeated on demand with different settings. CoGeBlast takes sequence from any other CoGe applications or text as query to any number of genomes; the blastn or tblastx results may be downloaded into GEvo panels. GEvo panels may be combined via links to create experiments.

Ks values may be calculated for each data point in SynMap if the genomes being compared are repeat-masked and have annotated CDS sequences. The 4a.53 maize sequence was used. Several genomic comparisons in SynMap have Ks values pre-calculated, including sorghum/maize. Syntenic gene pairs were identified by using blastn with SynMaps's default settings [−W (word size) = 11, −G (gap open) = −1, −E (gap extend) = −1, −q (mismatch) = −3, −r (match) = 1] and an e-value cutoff 0.05. These pairs were used to identify any putative homeologs between coding sequences using DAGchainer to identify collinear sets of putative genes with the following parameters: −D = 20, −g = 10, −A = 5. Ks values for syntenic gene pairs were calculated by first performing a global alignment of virtual protein sequences using the Needleman-Wunsch algorithm [31] implemented in python (http://python.org/pypi/nwalign/). The BLOSOM62 scoring matrix was used for the alignments [32]. From these protein alignments, the codon DNA alignment was generated through back-translation. Ks values were calculated using codeml of the PAML software package [33] on the codon alignment with the following parameters: outfile = mlc, aaDist = 0, verbose = 0, noisy = 0, RateAncestor = 1, kappa = 2, model = 0, ndata = 1, aaRatefile = wag.dat, Small_Diff = .5e-6, CodonFreq = 2, runmode = −2, alpha = 0, omega = 0.4, fix_kappa = 0, Mgene = 0, method = 0, fix_omega = 0, getSE = 0, NSsites = 0, seqtype = 1, cleandata = 0, icode = 0, fix_alpha = 1, clock = 0, ncatG = 1, Malpha = 0, fix_blength = 0. This pipeline is part of the SynMap application in the CoGe suite of comparative genomics software, and its dotplot visualization tool was used to generate the Ks color-coded lines of Dataset S2, and its text output was used to supply the Ks values for the “sorghum/maize Ks differences” Methods section to follow. Since this Ks pipeline will calculate Ks values for erroneously aligned pairs, values far off from an expected normal distribution for any experiment were discarded.

Preparing the Sb-Zm1-Zm2 Sequences for Comparison

The entire sorghum genome was subjected to a 50× repeat mask, where every base pair that was covered more than 50 times by a blast hit from a whole-genome self-self blast was masked, using parameters of blastn at word size 16 and e-value cutoff of 0.001. Repeats over 50× genome-wide were masked by changing their sequence to “x.” We needed to use a step-wise approach to accomplish the same 50× mask for maize because a direct self-blast was too memory-intensive for our computers. First we self-blasted pseudomolecules 1–3 as if they were the whole maize genome and masked their 50× repeats. Then we added these 3 larger masked chromosomes to the other 7 unmasked and performed self-BLAST—as with sorghum above. Repeated sequences are color-coded pink in Figure 3, panels B and C. The sorghum 50× masked genome was further masked for every sequence that is not either an Sbi1.4 exon or other sequence shared orthologously with rice, as derived from the rice/sorghum Pipeline 1.0. The non-exon, non-conserved sorghum sequences masked by this method are colored orange in our GEvo graphics (e.g. Figure 3, panel A).

The nine Sb-Zm1-Zm2 regions were derived from SynMap blastn [34] dotplots using the DAGchainer settings −g = 10 genes, −D = 20 genes, −A = 5 genes, and a Ks color code that clearly distinguishes syntenic lines reflecting sorghum/maize orthologs to lines reflecting more ancient syntenies. When a single stretch of sorghum clearly hits two longer stretches of maize, the center of the overlapped region was used as an anchor to create Sb-Zm1 and Sb-Zm2 GEvo panels, which are then combined into a single view. The sorghum/maize Ks-colorized dotplot can be seen in Dataset S2, where the identification of Sb1 is illustrated. It is possible to regenerate a near-identical graphic in CoGe by visiting http://tinyurl.com/ygx2apu. The 28 additional Sb-Zm1-Zm2 regions were discovered by choosing as query exons from sorghum genes that encode transcription factors. Each query found, using CoGeBlast, two orthologs in maize about one-third of the time. From CoGeBlast output, it is easy to create Sb-Zm1-Zm2 GEvo panels. Lengths of these three chromosomes were adjusted so that a chosen segment of sorghum begins and ends with a retained gene, was entirely represented syntenically within the two maize segments, and syntenic coverage did not improve by adding 500 kb on both sides of the maize chromosomes. Inversions do not cloud our analyses because all inversions we include begin and end within each region.

Our primary data of Dataset S1 required that every gene on the sorghum gene list receive one among several possible annotations. Genes in local arrays were marked as parent, duplicate (D or DUP), or interrupter (a gene located within a tandem repeat) using published methods [7] and duplicates were marked and ignored subsequently; up to three interrupter genes were permitted. If a remaining gene occurred syntenically (blastn bitscore >50) on a maize homeolog, then it was coded “1” or “2” if it occurred on only one of the homeologs or “B” if it occurred on both. A few genes were invalidated for technical reasons, and some genes were not found in the syntenic position in either maize homeolog (encoded as “0”). Genes represented by fragments were counted as “present” even though they were almost certainly in the process of removal. In this manner, each of our 37 Sb-Zm1-Zm2 regions were reduced to a code of shared genes, like B1122BBB12BB2B21BBB11B21B, and trimmed to begin and end with a B (present in both maize homeologs) where the terminal Bs were not within inversions. For the diagrams of Figure 3 B, C, and D and Dataset S2, and for all analyses of runs, as discussed in the text we removed runs of 1's or 2's that extended beyond nine genes. This is because our analyses suggested that a run of 10 or more 1's indicates that the 10 genes that would be the corresponding 2's had jumped elsewhere in the genome. The unmodified data are in Dataset S1. At this time, accurate fractionation annotations would be difficult or impossible to achieve automatically largely because of biological complications involving inversions and also by contig misalignments during sequence assembly.

Fractionation Statistics

The binomial test was used to evaluate the probability that the ratio of deletions on the maize homeologs could occur by chance given an expectation that a single deletion is equally likely to occur on either homeolog. The distribution of all observed deletion lengths is plotted in Figure 4 as the blue bars for the over- and under-fractionated homeologs. Using the initial hypothesis of a deletion mechanism that independently eliminates one gene at a time, a simulation of gene loss was carried out. Starting with a length equal to all genes, both deleted and still present, genes were deleted at random until the simulated number of deletions was equal to the true observed number. The distribution of apparent deletion lengths for the run was then saved and the preceding steps were repeated 1,000 times. This gives a distribution of frequencies of all deletion lengths. The median number of apparent deletion runs from these simulations is shown by the white circles in the grey lines of Figure 5, with grey line itself marking the values between which the results from 95% of the simulations fall.

For Figure 5E, which plots runs of genes conserved on both maize homeologs, the above model was modified by generating two lengths each equal to the total number of sorghum genes within the dataset, and then deleting genes from either one or the other sequence (with an bias for deleting genes from one or the other dataset equal to that observed in the overall fractionation dataset) until the number of retrain genes (Bs) was equal to the true number observed, with the constraint that once a simulated gene was deleted from one dataset, the orthologous gene in the other dataset would never be deleted.

As the simulated distribution did not perfectly match the observed results, a genetic algorithm using 20 (genetic) character states, each representing a 5% (1/20) chance that a deletion would be some length between 1 and 5 genes long was used to determine, given the region length and the distribution of observed deletion lengths, the ratio between different deletion lengths to use in the simulation described above to achieve the best match between simulated and observed data. The fitness of solutions in the evolutionary algorithm were scored using the Monte Carlo method described in the proceeding paragraph (with the modification that rather than fixing the deletion length at 1 gene, deletion lengths were selected using the weighted averages generated by the evolutionary algorithm) with the most fit solutions being those where the median simulated number of deletion runs was least different from the observed number of runs. The genetic algorithm was allowed to run for 100,000 generations. These new weighted average deletion lengths can then be used to generate new sets of expectations for data, as seen in Figure 4D. The script used to run the genetic algorithm is available at http://code.google.com/p/bpbio/source/browse/trunk/scripts/fractionation/fractionation_ga.py and in Dataset S7.

Finding Segmental Translocations

Sorghum genes with known orthologs in rice were blasted against the sorghum (JGI 1.4 gene models), rice (TIGR 6.0), and maize (4a.53 filtered gene set; maizegenome.org) datasets. We used the score of the best sorghum-rice alignment as a cut-off to avoid hits from genes that diverged before the rice sorghum split, and removed genes with more than one hit above that threshold in the sorghum-sorghum blast to avoid the inclusion of genes that duplicated in the sorghum-maize lineage since the divergence from rice. These criteria left us with a set of approximately 10,000 genes with a single hit in sorghum that had a greater bit score than any hit in rice, and one or more hits satisfying the same conditions in maize. 406 genes from this dataset overlapped with genes identified as retained (noted as “B” in Dataset S1) by manual annotators, and 771 overlapped with genes identified as fractionated. Stretches of 10 or more genes deleted from the same chromosome were identified on Dataset S1 and the missing region was identified by a discontinuity in the appropriate sorghum/maize dot plot. We built a string of exons that identified each gene in the deleted region and used it as query to the subject maize genome. The maize genome was 50× repeat-masked, as described, and blastn used settings of word size 7, and e-value <0.001. Hits were achieved in CoGeBlast and evaluated in GEvo. Any three of the expected genes, arranged syntenically, in unexpected regions of genome were taken as evidence for a segmental translocation even though a gene might have been represented by a fragment rather than an entire gene.

Evaluation of Copy Number of Genes in the Nuclear Genomes of Sorghum and Maize

The coding sequences of the subset of genes from the JGI sorghum 1.4 gene set that had been identified as orthologous to a single rice gene were blasted against the MSU6 rice gene set and the maize 4.a53 filtered gene sets as well as against the same sorghum gene set. For each sorghum gene, the bit score of the highest-scoring alignment against a rice gene was used as a cutoff to exclude hits from genes that had diverged from the gene being tested before the rice/sorghum split. Sorghum genes that hit one or more additional sorghum genes with bit scores higher than that cutoff were excluded from the analysis to exclude genes duplicated in the maize/sorghum lineage since the rice/sorghum split.

The number of hits to genes in the maize filtered gene set for the remaining sorghum genes (with scores higher than the best hit in rice) was recorded. After the accuracy of a sample of the results were manually checked using CoGe, the final data were generated by looking at the average number of maize genes found using this process for genes assigned to the fractionated and unfractionated categories by manual annotation.

Looking for Sb-Zm Ks Differences Depending on Whether Over- or Under-Fractionated

Ks values for shared open reading frames in sorghum and maize (4.a53) were precalculated and loaded into SynMap, in CoGe as described previously in METHODS. The sorghum/maize orthologs that also fell into the 37 regions that were hand-annotated for the primary fractionation data (Dataset S1) were identified. Next, sorghum genes that hit to genes in both maize homeologs (encoded “B”) were paired and their Ks values compared. Data were reported in the format Sb-Zm1-under-fractionated/Sb-Zm2-over-fractionated. Visual examination of the Ks data showed a minority portion of very extreme ratios, likely the result of misalignments, alignments to pseudogenes, or alignments to non-orthologous genes. Such misalignments were expected due to the fragmented nature of many B73 genes and contig assembly error. The 16% of pairs with the most extreme ratios as compared to the median were removed from the dataset and not used to calculate results.

Counting Methylation Domains on Maize Chromosomes

We overlaid McrBC methylation data from [22] onto the annotated maize pseudomolecule sequence (dataset Zm 4a.53) and uploaded the modified database into the genome viewer we use with CoGe: GenomeView. We were able to visualize on GenomeView the locations of methylated sites on maize chromosome regions. After anchoring both maize homeologs to their orthologous sorghum sequence with the stop-start sites used in our fractionation analyses, we manually counted the number of methylation peaks in each maize homeologous region in question.

Analyses of Deletions within Exons

Using GEvo, with parameters set for blastn with a spike-length of 15 bp, we visually scanned all retained maize genes from our Sb/Zm1/Zm2 dataset to look for gaps within exons of one or the other maize homeolog. This level of resolution did not permit us to identify single gaps less than approximately 15 bp long. However, we did not intend to be exhaustive. Once a gap was identified, we extracted the sorghum exon sequence and used it as query in a blastn comparison to rice; this use of the rice as a secondary outgroup often confirmed the sorghum full-length exon annotation, and when it did, we re-blastn'd this sequence against the multiple subjects rice (Oryza sativa v5 masked repeats 50×X), sorghum (vSbi1.4 exons, 50×X mask+syntenic thread with Os), and maize v4a.53 to produce GEvo images like that shown in Figure 5A. We then took the corresponding exon sequence data from rice, sorghum, and both maize homeologs and used ClustalW (http://www.ch.embnet.org/software/ClustalW.html) to visualize the sequence alignment surrounding the gap, as well as the sequence on the homeolog without the deletion (as in Figure 5B).