Replicated high-density genetic maps of two great tit populations reveal fine-scale genomic departures from sex-equal recombination rates.

Linking variation in quantitative traits to variation in the genome is an important, but challenging task in the study of life-history evolution. Linkage maps provide a valuable tool for the unravelling of such trait-gene associations. Moreover, they give insight into recombination landscapes and between-species karyotype evolution. Here we used genotype data, generated from a 10k single-nucleotide polymorphism (SNP) chip, of over 2000 individuals to produce high-density linkage maps of the great tit (Parus major), a passerine bird that serves as a model species for ecological and evolutionary questions. We created independent maps from two distinct populations: a captive F2-cross from The Netherlands (NL) and a wild population from the United Kingdom (UK). The two maps contained 6554 SNPs in 32 linkage groups, spanning 2010 cM and 1917 cM for the NL and UK populations, respectively, and were similar in size and marker order. Subtle levels of heterochiasmy within and between chromosomes were remarkably consistent between the populations, suggesting that the local departures from sex-equal recombination rates have evolved. This key and surprising result would have been impossible to detect if only one population was mapped. A comparison with zebra finch Taeniopygia guttata, chicken Gallus gallus and the green anole lizard Anolis carolinensis genomes provided further insight into the evolution of avian karyotypes.

Introduction

In the last decade, linkage mapping projects in the bird species have been extended from studies of agriculturally important organisms, such as chickens (Groenen ), to a much more diverse set of species including other galliformes such as the turkey (Reed ), anseriformes such as the mallard duck (Huang ) and several passerine species including the collared flycatcher (Backström ), the blue tit (Hansson ), the zebra finch (Stapley ; Backström ) and the great reed warbler (Hansson ). Although the most comprehensive whole-genome maps have been constructed for galliformes it is notable that passerines, by far the largest order in the class Aves, represent the majority of other species for which such resources have been developed. Linkage maps have been constructed with two main purposes in mind: (1) to identify regions of the genome that are responsible for genetic variation in quantitative traits, that is, to map quantitative trait loci and (2) to understand the genomic features and evolution of avian karyotypes. Often, however, particular species are well suited to investigating both types of question. In this manuscript we describe linkage maps of the great tit (Parus major) a species that is widely studied in ecology and evolutionary biology (see below).

A motivation for constructing linkage maps in wild bird species, especially passerines, which have been the focus of many ecological and evolutionary studies (Bennett and Owens, 2002), is to provide understanding of how evolution operates in natural populations. Quantitative trait loci mapping can unravel knowledge of the molecular genetic basis of ecologically relevant traits, which are usually expected to have a complex genetic architecture (Slate, 2005). Although linkage mapping is being superseded by genome-wide association studies, an approach that does not directly require a linkage map, linkage maps remain useful tools for identifying the chromosomal location and order of markers used in a genome-wide association studies.

In addition to understanding the genetic basis of ecologically important traits, linkage maps provide valuable insight into the evolution of genome organisation. For example, they help to reveal and explain the recombination landscape across the genome, how recombination rates might differ between the sexes and how species are characterised by structural genomic changes. Studies of bird genome evolution have identified a number of interesting features. First, bird karyotypes are remarkably highly conserved, with the majority of species having a haploid number of 38–42 chromosomes, and most species showing conserved synteny (Kayang ; Backström ; Stapley ; Ellegren, 2010). Second, recombination is usually not evenly distributed between and within chromosomes. Macrochromosomes tend to have lower recombination rates (per Mbp) than microchromosomes, presumably because all chromosomes have an obligate crossing-over during meiosis, even if they are physically short (Groenen ; Stapley ). Within macrochromosomes of zebra finches, recombination rates are lower at the centre relative to telomeric regions (Stapley ; Backström ; Stapley ), although this phenomenon has not been observed to the same extent in chickens (Groenen ). Third, some bird species display heterochiasmy, that is, a difference in recombination rates between sexes.

Heterochiasmy was first believed to reflect pleiotropic consequences of selection against recombination between the sex chromosomes (the Haldane−Huxley rule) (Haldane, 1922; Huxley, 1928). More recently, several evolutionary hypotheses have been developed, where generally differences in recombination rates between sexes are expected to be due to sex differences in selection pressures. For example, in mammals recombination in males is commonly found to be lower than in females (Popa ), possibly to reduce the possibility of recombination breaking up favourable gene combinations in the sex with more intense selection (Trivers, 1988). Other possible evolutionary explanations for heterochiasmy are that the sexes have different selective optima to balance the deleterious mutation rate and the cost of efficient repair (Ellegren, 2007). However, pronounced heterochiasmy is not a general feature in all animals and is certainly not ubiquitous in birds (Stapley ; Groenen ). The most striking example has been observed in the great reed warbler, where a female:male map length ratio of 2.15 was initially observed (Hansson ), although a subsequent study with more markers revised this down (female:male ratio=1.56; Åkesson ). In the blue tit recombination rates are slightly greater in females (female:male ratio=1.07; Hansson ) and in a study of turkey the opposite pattern was found: (female:male ratio=0.86; Aslam ). Apart from the overall differences in total genome map length for males and females, substantial variation in recombination rates between sexes may also exist within and between autosomes of the same species (Dunn and Bennett, 1967; Robertson, 1984). For example, in humans, levels of heterochiasmy are associated with the physical location on the chromosome relative to the centromere, where recombination rates tended to be female biased near the centromere and male biased towards the telomeres (Broman ), with a possible relationship between heterochiasmy and GC content (Popa ). In birds, levels of heterochiasmy have not been studied at either the within- or between-chromosome scale, which are essential starting points if trying to understand how sexual selection may cause variation in recombination rates (Ellegren, 2007).

In most species, especially of natural populations, mapping studies have been conducted on single populations only (but see Poissant for an exception, albeit with a modest number of markers). In order to make further progress in understanding bird genomes, high-density linkage maps, replicated across populations of the same species are highly desirable. The initial observation that gene order was highly conserved over large evolutionary distances, for example, between passerines and galliformes, was shown to be only partially true once higher density maps with sufficient resolution to detect minor within-chromosome rearrangements became available (Backström ; Stapley ; Backström ). However, alternative marker orders often have similar statistical probabilities (likelihoods) and so determining whether subtle differences in gene order are real is not straightforward. If a particular marker order is identified in two independent populations of the same species, the evidence for rearrangements becomes much greater. Similarly, the replication of linkage maps in two populations can help understand whether heterochiasmy is present on a (genomically) local scale. If particular regions of the genome display apparent heterochiasmy, then without replication, it is difficult to rule out random, non-repeatable, sex-biased sampling of the location of recombination events as the cause, rather than a conserved feature of a species' genome. Constructing maps in two independent populations overcome this problem, because observing similar patterns in both populations provide strong support for true species-level differences in male and female recombination rates.

The great tit, a passerine of the family Paridae, is a model organism for evolutionary and ecological research in wild populations, with a wide body of work, addressing topics such as life history evolution (Garant ; Postma and van Noordwijk, 2005), parasite infestation (Norris ; Heeb ), climate change (Visser ; Charmantier ), plumage characteristics (Evans and Sheldon, 2012), song (Slabbekoorn and Peet, 2003; Naguib ) and personality (Dingemanse ; Drent ; Quinn ; van Oers and Naguib, 2013). The development of a genome-wide linkage map is, therefore, essential to help shed light on the link between genotypic and phenotypic variation in such studies. Using chromosome painting, the karyotype of great tit has been estimated at 2n=80 (Nanda ), similar to other passerines including the zebra finch (Itoh and Arnold, 2005). The zebra finch and great tit diverged from a common ancestor 40–45 million years ago (Barker ), but given the highly conserved nature of avian karyotypes, a high degree of shared synteny between zebra finch and great tit genomes is expected. As with most bird species, great tit chromosomes vary in size (Nanda ), over a range of almost two orders of magnitude. Around 10 large macrochromosomes typically have similar sizes to most mammal chromosomes and in addition there are numerous smaller microchromosomes (Ellegren, 2010).

In this study, we describe the construction and analysis of two high-density whole-genome linkage maps of the great tit. To construct the maps we used a set of SNP markers originating from a 9193 SNP chip (van Bers ). We genotyped a captive F2 inter-cross population from selected lines on high and low early exploratory behaviour (Drent ; van Oers ) at the Netherlands Institute of Ecology (NIOO-KNAW), NL and birds from a wild pedigreed long-term study population (Wytham Woods) near Oxford, UK. Producing independent replicate maps for separate populations should (i) give insight into how accurate/robust linkage maps are, (ii) allow tests of whether variation in recombination rates across a chromosome are consistent across populations, (iii) allow tests of whether small differences in recombination rates between sexes (heterochiasmy) are chance departures from 1:1 or are consistent across populations and (iv) allow tests of whether minor rearrangements between bird species are genuine or an artefact of minor errors caused by alternative marker orders having very similar likelihoods. It will, thus, provide more robust insight into avian karyotype evolution, if within-species patterns are demonstrated to be repeatable.

Materials and methods

SNP genotyping

A total of 4702 birds were genotyped with an Illumina iSelect BeadChip containing 9193 SNPs, according to the manufacturers' protocol (van Bers ). Genotyping results were analysed using Genome Studio software (Illumina, San Diego, CA, USA). A total of 8340 SNPs were successfully genotyped in both the NL and the UK population, but only 7215 and 7203 of these SNPs were polymorphic in the NL and UK population, respectively. After excluding markers with very low minor allele frequency (MAF<0.01), 6897 and 7032 SNPs were used for map construction in the NL and UK populations, respectively. Detailed descriptions of SNP identification, SNP selection and genotyping results are given elsewhere (van Bers ; Santure ; van Bers ).

Study species and mapping populations

For the construction of the linkage maps, we used 398 great tits originating from a captive population at the NIOO-KNAW, Wageningen, NL and 1656 individuals from a pedigreed long-term natural population (Wytham Woods) near Oxford, UK. The genotyped individuals in the NL population were selected from 480 genotyped individuals of the fourth generation of artificial selection on fast and slow early exploratory behaviour (Drent ) and the F1 and F2 intercrosses between the fast and slow lines (van Oers ). The final mapping pedigree (N=398) included all F2 intercross individuals (N=251), their F1 parents, the siblings of these F1 parents and all full-sib selection line families with their parents. These individuals were manually split into 35 subfamilies, with mean family size 15.0 (range 5–35). Sex and parentage of all individuals was confirmed from the genotypes. Owing to a small number of parent-offspring genotype mismatches, we removed five individuals from the data set, leaving 393 individuals for analysis (Supplementary Methods).

A subset of 1656 Wytham individuals with confirmed first degree pedigree links (genotyped parents and/or sibs and/or offspring) and confirmed sex were chosen for constructing the UK great tit linkage map (see Supplementary Methods for more details). For the UK population, complex pedigrees were first split into 61 subfamilies (of up to approximately 50 individuals and three generations) using the CRIGEN command (options: -size 50 and -gen 3) in Crimap (see below). Parent-offspring genotype mismatches were removed from the data set.

Predicted locations by comparative mapping

To identify and understand rearrangements in avian genomes, it is necessary to compare the positions of the homologues of typed markers on different species' genomes. Recently, the genome of the green anole lizard A. carolinensis was sequenced and assembled (Alfoldi ). Here we examined chromosomal rearrangements between great tit, zebra finch, chicken and the green anole; the latter provides a useful out-group to the avian species. The position of the great tit SNPs on the zebra finch, chicken and green anole lizard genomes was predicted based on the alignment of 23 to 150 bases of flanking sequence to the respective genome sequence using the default parameters in blastn (with setting a minimum expectation value (e-value) of 10−5). The three genome sequences used for comparative mapping were downloaded from the Ensembl Genome Browser; zebra finch: version Taeniopygia_guttata.taeGut3.2.4.61, www.ensembl.org/Taeniopygia_guttata/Info/Index; chicken: version Gallus_gallus.WASHUC2.56, www.ensembl.org/Gallus_gallus/Info/Index and green anole lizard: version Anolis_carolinensis.AnoCar2.0.61, www.ensembl.org/Anolis_carolinensis/Info/Index. Predicted physical map locations are the midpoints between the start and end of the mapping alignment and are, therefore, expected to be within 100 bases of the actual SNP position. Linkage groups and putative chromosomes were orientated according to the orientation of the chicken chromosomes.

Map construction

For linkage mapping we used a modified version (2.503) of the CriMap software (Green ). We chose to present maps with three different marker densities in this manuscript: (i) framework maps, incorporating markers that could be assigned a most likely order with LOD>3, that is, the best order is 103=1000 times more likely than the next best alternative order (ii) comprehensive maps, where ‘non-framework' markers were added by lowering levels of stringency (LOD>0.1) and (iii) parsimonious maps, adding as many markers as possible to the comprehensive map with LOD<0.1 (Supplementary Methods).

We applied two different initial strategies for linkage group construction. For the NL map, we assigned 6877 SNPs to putative great tit chromosomes based on sequence homology between SNP flanking sequences and the zebra finch genome sequence. For the UK map, we used the TWOPOINT output of CriMap to create linkage groups with the AUTOGROUP command (see Supplementary Methods for more details). In both mapping populations, the second allele of SNPs on the sex chromosome (Z), was assigned a dummy value (6) in females, to account for the hemizygosity of the female Z chromosome.

Heterochiasmy

Heterochiasmy was examined on three different scales, using the sex-specific maps option in CriMap (see section Map Construction). First, the genome-wide male and female total map lengths were compared. Second, male and female map lengths were estimated for each chromosome. Third, local levels of heterochiasmy were estimated using sex-specific map lengths in a window size of 20 markers and then sliding the window in five marker increments. For this analysis, markers common to comprehensive maps in both populations were used: Parus major (PMA) linkage groups PMA1, 1A, 2–4, 4A, 5–15, 17–21, 23, 24 and 26 contained enough shared comprehensive map markers to estimate local heterochiasmy. Following Poissant , heterochiasmy was measured using the size dimorphism index (SDI) developed by Lovich and Gibbons (1992). This index is considered the best estimator of sexual size dimorphism because, in contrast to ratios, it does not suffer from the fact that heterochiasmy is bound at zero in one sex and unbounded in the other sex. The SDI was obtained with the following formula:

Following convention (Lovich and Gibbons, 1992; Poissant ), estimates were then made positive when the female value was largest and negative when the male value was largest. To distinguish between whole-chromosome/whole-genome heterochiasmy and local heterochiasmy, we use the notation SDI and SDI_local, respectively.

Statistical analysis

We used Pearson's correlations to compare chromosome sizes between populations. We tested for the presence of genome-wide sexual dimorphism using sign tests, with the SDI of each chromosome being one data point, and tested the relationship in SDI between populations using reduced major axis regression (NLRA in PASW Statistics for Mac Version 18.0.3, SPPS Inc., Chicago, IL, USA). To calculate whether linkage group and interval SDI values fell outside the distribution, we conducted an outlier analysis. The critical values to treat a data point as an outlier were calculated by multiplying the difference between the 25% (Q1) and the 75% (Q3) quartiles by 2.2 and adding this to Q3 or subtracting it from Q1 (Hoaglin ). Data points were classified as outliers, when they were either larger than the Q3 or small than the Q1 critical value. To analyse whether local heterochiasmy was conserved between the two populations the auto- and cross-correlation functions implemented in R version 2.13 (R Development Core Team, 2011) were used. All other statistical analyses were conducted using PASW Statistics for Mac version 18.0.3 (SPSS Inc.).

Results

Initial assignment of markers into linkage groups

For the NL population (N=393 birds), 6877 markers (mean informative meiosis±s.d.=224.6±94.4) were grouped into 32 putative great tit chromosomes based on the alignment of markers with the zebra finch genome. In the UK population (N=1656 birds), 6094 markers (mean informative meiosis±s.d.=408.0±124.2) were grouped into the same 32 putative great tit chromosomes, using AUTOGROUP on the set of autosomal SNPs. The grouping of markers by the AUTOGROUP option in CriMap demonstrated a remarkable conservation of synteny between the great tit and zebra finch genomes. Based on the very high proportion of markers within each group that aligned to a single zebra finch chromosome (range 91–100% mean=99.7), chromosome identities could be assigned unambiguously to each linkage group from homology with the zebra finch genome. We will, therefore, refer to the homologous chromosome number when referring to the linkage groups throughout the paper. This follows the convention of Stapley , who assigned zebra finch linkage group names based on homology with the chicken assembly. Thus, great tit, zebra finch and chicken chromosome names are consistent. This method is in contrast to the historic approach of naming chromosomes in order of decreasing physical size. For example, from fluorescence in situ hybridization paints in great tits, PMA2 is known to be larger than PMA1 (Nanda ).

Linkage map characteristics

The NL and UK framework maps (LOD⩾3) contained 1766 and 1674 markers, respectively. Of these, 898 markers were present on both maps. When comparing the two maps, 37 out of these 898 markers had a position on both maps that could not be predicted from the position from the other map (Supplementary Figures 1 and 2). The mean interval between the position of these markers on the NL map and the position based on the UK map was 1.91 cM (range 0.02–5.40 cM).

The NL and UK framework maps spanned 1869 cM and 1893 cM, respectively, across 32 linkage groups, homologous to TGU1 to TGU15, TGU17 to TGU28, TGUZ and the zebra finch linkage group LGE22 (Table 1). These linkage groups are homologous to GGA1 to GGA15, GGA17 to GGA28, GGAZ and two chicken linkage groups, LG2 and LGE22C19W28_E50C23. Markers that aligned with chicken chromosomes GGA1, GGA4 and GGA25 were assembled into two distinct linkage groups (PMA1 and PMA1A, PMA4 and PMA4A, and PMA25A and PMA25B) in both great tit populations. Merging the paired linkage groups into one resulted in a distance larger than 50 cM between the two groups of markers. For PMA25A/PMA25B, the two markers that were positioned closest together on the zebra finch genome, but were assigned to different linkage groups (PMA25A and PMA25B) were only 150k base pairs apart on the zebra finch genome. One pair of markers, M3554 from PMA25A and M6440 from PMA25B did have a twopoint LOD score of 2.2 in the NL population, but there was no evidence for linkage between this marker pair in the UK population. Furthermore, no other markers pointed to a possible linkage between PMA25A and PMA25B, suggesting that the two linkage groups may represent two different chromosomes.

Table 1

Characteristics of the parsimonious great tit linkage groups (PMA) mapped in the NL and UK populations

PMA	NL				UK
	No. of markers	Map length (cM)	Female: male ratio	% coverage	No. of markers	Map length (cM)	Female: male ratio	% coverage
1	566	138.15	0.86	99.0	523	139.88	0.90	99.0
1A	433	112.82	0.96	97.4	334	93.58	0.91	97.8
2	664	127.41	0.60	93.0	624	139.69	0.72	99.1
3	602	107.05	0.65	95.3	525	114.87	0.73	99.8
4	370	97.33	0.78	99.2	306	97.65	0.90	99.2
4A	115	65.52	1.19	99.3	100	59.38	1.33	99.3
5	371	101.05	0.92	99.8	280	98.57	0.94	98.2
6	219	78.07	1.10	98.2	141	78.02	1.17	98.8
7	229	83.67	1.05	99.8	166	72.64	0.99	99.8
8	151	49.74	1.03	98.7	125	53.83	0.94	98.7
9	135	56.57	1.06	97.3	121	54.20	1.18	97.2
10	152	56.67	1.09	98.4	129	50.47	1.08	97.8
11	151	58.49	0.88	98.2	112	58.23	1.08	98.2
12	178	55.42	0.95	99.0	135	51.92	0.98	99.0
13	140	54.17	0.65	97.9	110	40.95	0.83	96.5
14	145	51.65	1.01	97.0	113	49.20	1.17	97.0
15	205	48.93	1.23	96.4	162	49.15	1.14	98.1
17	107	50.05	1.08	90.2	90	45.36	1.00	87.7
18	105	51.02	0.93	98.7	87	49.91	0.94	98.5
19	101	51.52	1.01	95.1	94	49.43	1.02	95.5
20	185	50.94	1.15	99.6	137	49.44	1.05	98.1
21	56	53.08	0.95	91.9	51	46.76	1.00	91.9
22	24	55.75	1.05	91.5	17	57.69	0.97	91.5
23	47	40.33	1.15	97.3	45	51.43	1.10	75.4
24	65	52.33	0.82	97.0	51	49.51	1.21	96.2
25A	5	4.58	1.48	96.6a	3	4.43	1.29	99.2a
25B	8	5.76	0.04	—	7	3.68	0.06	—
26	55	54.91	1.01	97.2	38	53.45	0.86	97.2
27	37	54.32	0.70	83.3	30	49.15	0.87	89.7
28	50	52.43	1.55	96.6	42	48.70	1.14	96.6
Z	179	55.53	NA	99.5	177	51.21	NA	99.5
LGE22	5	34.59	1.55	77.7	3	4.44	1.00	28.4
Total	5855	2009.85	1.07	97.1b	4878	1916.82	1.04	98.3b

Abbreviations: cM, centiMorgan; NL, Netherland; UK, United Kingdom.

The number of markers for each linkage group, the length (cM) of the sex-averaged maps, ratio of female-to-male map length of the sex-specific maps and the percentage of genome coverage are all reported. Genome coverage was calculated as the percentage of each zebra finch chromosome covered by the two most distal SNPs on the map.

For the calculation of the genome coverage of PMA25A and PMA25B, we combined the two linkage groups.

Total genome coverage is calculated over the whole genome and not as the mean of the separate linkage groups.

The comprehensive maps contained 2853 markers spanning 1867 cM and 3162 markers spanning 1889 cM for the NL and UK populations, respectively. For the parsimonious map, a total of 5855 SNPs were placed on the NL map spanning 2010 cM (Supplementary Figure 1), whereas the UK map consisted of 4878 SNPs spanning 1917 cM (Supplementary Figure 2). This corresponds to an increase in size of 7.7% (NL) and 1.5% (UK) compared with the framework maps. In total 6554 unique SNP markers were present on either the NL map or the UK map. Of these, 4189 were present on both the NL and the UK map, whereas 1666 SNPs were unique to the NL map and 699 markers were present on the UK map and not on the NL map. The MAF for the 1666 markers present on the NL map, but not on the UK map (NL specific: mean±s.e.m.; 0.339±0.003), was significantly lower (t5140=4.16, P<0.001) compared with the MAF of markers present on both maps (present on both: mean±s.e.m.; 0.352±0.001). The MAF for the markers that were present on the UK map, but not on NL map did not have significantly lower MAF (UK specific: mean±s.e.m.; 0.345±0.003) than those on both maps. The number of SNPs per linkage group ranged from 5 (PMALGE22) to 664 (PMA2) for the NL population and from 3 (PMALGE22) to 624 (PMA2) for the UK population. There were no marker inconsistencies at linkage group level (that is, no marker was present on a specific linkage group on one map and present on another linkage group on the other map). Within the linkage groups, we found high consistency between the marker orders in the two populations (Supplementary Figure 3), the only rearrangements being apparent were local inversions involving a small number of markers (see discussion). Recombination rates (measured as cM/Mbp) tended to be greater at the ends of macrochromosomes than in the central parts, although unlike the situation in captive zebra finch populations, recombination ‘deserts' were not observed in the middle of macrochromosomes (Supplementary Figure 4).

Linkage group sizes ranged from 4.58 cM to 138.15 cM in the NL population and from 3.68 cM to 139.88 cM in the UK population (Table 1), and these sizes were highly correlated between the populations (Figure 1; r=0.97, n=32; P<0.001), indicating that the recombination rates were not population dependent, and did not differ systematically between the wild and recently-founded captive populations studied here. Based on the alignment of the flanking sequences of the most distal SNPs in great tit linkage map to the zebra finch genome, the total genome coverage was 97.1% (range 77.7–99.8 for individual linkage groups) for the NL map and 98.3% (range 28.4–99.8) for the UK map (Table 1). Assuming that our linkage map covers 32 of the expected 40 chromosomes (that is, taking PMA25A, PMA25B and LGE22 as a separate great tit chromosome), and that these three and the eight missing microchromosomes have at least one obligate recombination event during meiosis, implying a map size of at least 50 cM per chromosome, the total predicted map sizes are ∼2700 cM for the NL and ∼2600 cM for the UK populations.

Figure 1

Total map size for each chromosome of the NL maps plotted against the map length of the UK maps. The line shows the expected value if map lengths from both populations were identical.

Genome wide, there was no evidence for large sex differences in map length. We found no significant levels of heterochiasmy (SDI in the NL or the UK population (Table 1; NL: x̄=−0.070, 95% C.I.=−0.166−0.0069, t27=−1.61, P=0.12; UK: x̄=−0.007, 95% C.I.=−0.067−0.051, t27=−0.21, P=0.84). At the chromosomal level SDI levels correlated strongly between populations (reduced major axis regression, y=0.72x+0.04, r2=0.57, F1,26=37.07, P<0.001; Figure 2). An outlier analysis showed that the sex-specific linkage maps of three linkage groups in the NL population significantly differed from each other: males in PMA2, PMA3 and PMA13 had longer map lengths compared with females in the NL population. Although not significant, these three linkage groups also had the lowest SDI values (that is, greater recombination rates in males) in the UK population.

Figure 2

Comparison of levels of heterochiasmy per linkage group for the NL and the UK populations. Heterochiasmy is indicated as the size dimorphism index (SDI). Negative SDI values indicate larger male linkage groups compared with female linkage groups. The solid line represents the y=x line, the dotted line the outcome of the reduced major axis regression (y=0.72x+0.04). Only linkage groups larger than 10 cM are plotted. Linkage groups that deviate significantly from 1:1 are indicated with *.

Because replicated maps from independent populations were constructed, it was possible to examine whether subtle, more localised differences in heterochiasmy were conserved. Most revealingly, the plots (Figure 3) and analyses of SDI_local revealed fine-scale heterochiasmy that was remarkably conserved between the two populations. Cross-correlations, in SDI_local, between the two populations were positive at all 22 chromosomes considered (range in r=0.05–0.999; median r=0.668; Supplementary Table 1). A sign test on all 22 correlations being positive is highly significant (P=4.8 × 10−7, assuming a null expectation of a correlation being equally likely to be positive or negative). At 17 of the 22 chromosomes the cross-correlation would be regarded as significantly positive, if each window was regarded as an independent data point. Importantly, the overall pattern of local heterochiasmy being highly conserved across the two populations was not dependent on window size or the interval from one marker to the next (See Supplementary Figures 5a–c, for SDI_local plots using alternative window sizes).

Figure 3

Size Dimorphism Index (SDI) calculated for windows of 20 SNP markers on the NL (orange) and UK (blue) framework maps. Subsequent windows were chosen by sliding the window five SNP markers along the linkage group. Positive SDI indicates that female recombination rates within the 20 SNP window was higher; negative values indicate greater recombination rates in males. Only linkage groups with at least three windows are plotted.

Comparative mapping

For 5396 of the 5855 SNPs that mapped on both the NL and the UK parsimonious linkage map, we found homology to the zebra finch genome; 1231 of the SNPs could be assigned locations on the chicken genome and 169 were aligned to the green anole lizard genome.

In general, gene order was conserved within syntenic blocks when comparing the map positions on the great tit parsimonious map with the physical position on the zebra finch (Supplementary Figure 6) and the chicken genome (Supplementary Figure 7). Two linkage groups (PMA8 and PMA19) did not show any intrachromosomal rearrangements within either the zebra finch–great tit or the chicken–great tit comparison, that is, marker order was exactly the same in great tit, zebra finch and chicken (see example in Figure 4a). For four linkage groups (PMA5, PMA9, PMA12 and PMA14) chromosomal rearrangements were similar whether comparing great tit to zebra finch or to chicken (see Figure 4b). These cases represent rearrangements that have occurred in the great tit lineage since great tits and zebra finches diverged. In four cases, rearrangements were only visible for the TGU-PMA comparison (PMA2, PMA10, PMA17 and PMA18) but not the GGA-PMA comparison (see Figure 4c). These represent rearrangements in the lineage leading to zebra finches after the great tit-zebra finch divergence. In five linkage groups (PMA4, PMA6, PMA13, PMA15 and PMA20) rearrangements were present only in the GGA-PMA comparison (see Figure 4d), and represent rearrangements that occurred after galliformes and passerines diverged, but before the divergence of great tits and zebra finches. The remaining linkage groups have undergone inversions in the great tit compared with either chicken or zebra finch, but more complex rearrangements were also observable. For example, in PMA1, represented by TGU1 and the distal part of GGA1, two inversions and a third inversion within the second inversion were present in the TGU-PMA comparison, although a rearrangement of ∼30 Mb is also present in the GGA-PMA comparison (Figure 4e). In another example, PMA7, an inversion of ∼8 Mb is present in the GGA-PMA comparison, although in the TGU-PMA comparison a large inversion of ∼25 Mb is also present in the central region of the chromosome (Figure 4f). Inference of ancestral states is more complicated for these cases, even when high-density maps from multiple passerine species are available.

Figure 4

Comparative maps for six example great tit linkage groups, with the great tit linkage map positions plotted against the predicted physical position on the zebra finch (TGU; left panels) and chicken (GGA; right panels) genomes. (a) PMA19 is highly conserved across avian lineages; no rearrangements relative to TGU or GGA, (b) PMA14 is an example of a rearrangement specific to the great tit lineage; an inversion relative to marker order in both TGU and GGA, (c) PMA10 is an example of a TGU-specific rearrangement; an inversion between PMA and TGU, (d) PMA6 is an example of a difference between passerines and galliformes; rearrangement relative to GGA but not TGU, (e) and (f) are chromosomes that are generally evolutionarily less stable across avian genomes, where (e) PMA1; inversions compared with TGU, and both inversions and rearrangements compared with GGA, and (f) PMA7; distinct inversions when PMA is compared with TGU and GGA.

Out of the 32 great tit linkage groups, 22 contained at least one SNP that could be assigned a location on the green anole lizard genome. Thirteen of these linkage groups had shared syntenic segments, with at least two SNP markers linked on both a green anole lizard chromosome and a great tit linkage group (Figure 5). PMA1 and PMA1A aligned to different green anole lizard chromosomes, suggesting that a single chromosome 1/1A in chickens is a derived state.

Figure 5

Great tit-green anole lizard (Anolis carolinensis; ACA) synteny map. The horizontal bars represent the six A. carolinensis macrochromosomes ACA1–ACA6 and three microchromosomes LGa, LGb and LGc, to which aligned SNPs or blocks of SNPs could be ordered. Each coloured block represents a great tit (PMA) chromosome as is indicated in the lower part of the graph. Adjusted from (Alfoldi ).

Discussion

In this paper we show that genetic maps are highly conserved between two great tit populations, using the same genetic markers, but different types of mapping population and different map-building methods. Owing to the construction of replicated independent genetic maps we were able to detect differences in recombination rate between sexes at a scale that has not previously been described in any bird species. Independent linkage maps also meant that genuine chromosomal rearrangements between species could be distinguished from alternative marker orders with similar likelihoods.

According to the Haldane–Huxley rule, the recombination rate should be reduced in the heterogametic sex, that is, female maps are expected to be shorter than male maps in birds (Haldane, 1922; Huxley, 1928). Our results add to a number of birds where whole-genome map sizes between males and females do not differ, at least not dramatically (Stapley ; Groenen ). However, despite the relatively little genome-wide difference in male and female recombination rates, our results show that heterochiasmy is present at a more subtle, yet consistent scale both between and, especially, within chromosomes. Importantly, these findings are repeatable across populations, indicating that they are real features of the great tit genome rather than sex-biased random sampling of where recombination events were observed in the genome.

It is possible that the relative position of the centromere is important with respect to local heterochiasmy. Unfortunately, we do not have information on centromere positions, and can, therefore, only speculate on this. Several studies have shown that recombination rates in birds are lowest around the centromeres, and that they increase towards the telomeres (see for example, Backström ), though these patterns have not previously been shown to differ between the two sexes. One intriguing possibility that requires further exploration is that local variation in heterochiasmy is related to the location of sexually antagonistic loci. It is well known that directional selection on any trait can result in an increase in recombination rate (see Otto and Lenormand, 2002). Therefore, conserved heterochiasmy on a local scale could indicate regions of the genome where selection is stronger in one sex than the other. For example, loci that explain variation in male fitness might be expected to be in regions of male-biased recombination rate and vice-versa. Measuring the fitness consequences of individual loci in wild populations remains notoriously difficult, but is perhaps possible in populations with extensive life history data and rapidly developing genomics resources, such as the great tit. The discovery of subtle, yet real, sex differences in recombination rates illustrates the advantages of replicating mapping studies in two independent mapping populations.

Population differences

Overall, our results show that the linkage maps of the two populations are very similar. We found no evidence that positions of groups of consecutive markers differed between the two populations, indicating no gene order differences between the populations. We did find order differences between single markers in 10 linkage groups, which are most likely due to power issues in assigning correct orders on a small scale rather than true rearrangements between the populations. Slightly different marker orders typically produce very similar likelihoods and total map lengths when there are insufficient recombination events to identify the true order. In domestic chicken, Groenen compared maps between three mapping populations typed with the same SNP set. Though small differences were found, the authors explained this by differences in the data structure between the populations, causing variation in information content of the different markers, rather than assuming real population differences in gene order.

We used two different methods of assigning SNPs to chromosomes. AUTOGROUP and in silico mapping by BLASTing to the zebra finch genome produced near-identical linkage maps for the two populations, indicating that for avian species both these methods are robust, given a reasonably closely related genome for comparative mapping. Studies that cannot rely on reference genomes will, therefore, not be hampered by a lack of a priori information on linkage information when building maps and will most likely perform better when using the AUTOGROUP method. A disadvantage of using a reference genome for assigning markers to expected chromosomes is that blast errors may complicate the building process. Moreover, markers with low information content are unlikely to be included when assigning markers to linkage groups using AUTOGROUP, increasing build speed and accuracy. Overall, we were able to incorporate marginally more markers in the NL map (using BLAST) compared with the UK map (using AUTOGROUP).

Map lengths

The extrapolated complete genetic map sizes of 2600 and 2700 cM, are similar to map lengths of other passerines like the zebra finch (Stapley ; Backström ), blue tit Cyanistes caeruleus (Hansson ) and Siberian Jay Perisoreus infaustus (Jaari ), but they are smaller than those of chicken or turkey (∼3200 cM; Groenen ; Aslam ) and larger than map sizes of the great reed warbler (∼2000 cM; Åkesson ), though this map was of lower density, probably underestimating total map length. This further demonstrates that passerines have a reduced recombination relative to domesticated galliformes (Stapley ), as genome sizes of passerines are in general comparable with chicken, at around 1.1–1.2 Gbp (Bachmann ; Andrews ). Recombination is often not evenly distributed within chromosomes and in birds is often more pronounced near the telomeres. Therefore, map sizes could be underestimated when marker density is low and distal parts of chromosomes are poorly covered. The large number of markers used here (relative to other passerine maps), means that differences in map length between passerines and chickens cannot be attributed to differences in genome coverage. One possible explanation for differences in map length is that artificial selection during domestication has led to an increased recombination rate in chickens (Otto and Lenormand, 2002; Groenen ). We found slightly larger map lengths in our captive F2-cross NL population compared with the wild UK population; these differences were, however, small compared with the between species differences, and were possibly due to the higher number of markers in the NL map, leading to inclusion of lower quality markers rather than an effect of the few generations of selection. Previous studies have demonstrated that adding error-prone markers tends to increase overall genetic map lengths (Groenen , 2009; Ball ). Without maps of wild galliformes being constructed, it remains unknown whether the passerine-galliform difference is a result of domestication of chickens and turkeys or an ancestral difference between passerines and galliformes that predates domestication.

Interchromosomal rearrangements

We found only one putative interchromosomal rearrangement when comparing our map with other whole-genome high-density avian linkage maps. No previous evidence exists to expect that PMA25A and PMA25B represent two separate great tit chromosomes. A lack of power due to a lower number of markers on these linkage groups in combination with higher recombination rates in microchromosomes may, therefore, be the cause for not detecting linkage between these groups of markers. This observation is consistent with an earlier mapping study in zebra finch, where the homologues of several chicken microchromosomes were found to be split among more than one zebra finch linkage group (Stapley ), despite chromosome painting results, indicating that these linkage groups are physically on the same chromosome (Itoh and Arnold, 2005; Nanda ). However, when markers from PMA25A and PMA25B were mapped together in a single linkage group, there was relatively weak evidence of linkage between markers on these linkage groups in the NL data set and no evidence in the UK data set. At this stage it is unclear whether these linkage groups represent one or two chromosomes. Comparative mapping of the macrochromosomes gave identical results to a recent mapping study of the zebra finch, where markers aligning to GGA1 were found on three linkage groups (TGU1, TGU1A and TGU1B) and GGA4 was found to be homologous to two macrochromosomes (Stapley ). When more markers were added, TGU1 and TGU1B were merged (Backström ), which was confirmed in our study. This also confirms the results from chromosome painting studies for great tit, zebra finch and several other passerines (Itoh and Arnold, 2005; Nanda ).

The large blocks of synteny apparent by comparing the green anole lizard and the passerine great tit genome are, in general, consistent with an earlier comparison of green anole lizard and chicken (Alfoldi ). However, by comparing SNP locations in the great tit with their locations in the zebra finch, chicken and green anole lizard genomes, it was possible to infer the ancestral state of some of the largest differences between galliform and passerine karyotypes. Our results strongly suggest that chromosomes 1 and 1A were separate in a common avian ancestor, but became fused in the lineage that led to chickens. Similarly, PMA4 aligned to only one single Anolis carolinensis chromosome (A. carolinensis 5), whereas none of the 115 SNPs that were present on PMA4A aligned to a position on the A. carolinensis genome. Although this is not direct proof that PMA4 and PMA4A originate from two distinct ancestral chromosomes, it is consistent with the idea that these chromosomes also underwent a fusion in the chicken ancestor (Griffin ). The ancestral states of other chromosomes cannot be inferred from our data.

Intrachromosomal rearrangements

Despite the highly conserved synteny between chickens and passerines, we found numerous intrachromosomal rearrangements. The two replicated maps were very consistent in support for the rearrangements, allowing us to conclude that they are true species-level rearrangements. Inversions and other changes in gene order are considered to be an important factor in speciation processes as they can promote reproductive isolation (Noor ). A comparison of multiple galliform and passerine species is needed to shed further light on this aspect of avian chromosome evolution; this is an exciting area for further research, as we anticipate increasingly dense avian linkage maps to become available (see for example, Stapley ).

In conclusion, we have constructed high-density linkage maps of two independent great tit populations. The maps will be valuable resources to aid with quantitative trait loci mapping, genome-wide association studies and chromosome partitioning of quantitative genetic variation (Santure ; Robinson ). More fundamentally, we show that replicated maps constructed using the same SNPs reveal new insight into bird karyotype evolution; in particular a hitherto undetected degree of fine scale heterochiasmy that opens up an exciting new opportunity to study the evolution and fitness consequences of heterochiasmy.

Data archiving

Genotype and pedigree data along with SNP information were deposited in the Dryad repository: doi:10.5061/dryad.j7260.