Evidence for mostly non-parallel differentially methylated regions between urban and forest environments
We identified a total of 224 distinct DMRs between urban and forest great tits: 80 for Barcelona, 68 for Montpellier and 93 for Warsaw. Only 14 DMRs (6.25%) were found repeatedly in at least two comparisons, and only 3 were common to the three cities. 7 of these 14 parallel DMRs were in the same direction of methylation in urban compared to forest areas (Figure 5, Figure S4A). Barcelona urban birds presented significantly more hypomethylated than hypermethylated DMRs (χ2 = 11.25, P < 0.001), while no difference was found for Montpellier (χ2 = 0.941, P = 0.332) nor for Warsaw (χ2 = 0.011, P = 0.917). DMRs were distributed across all the 32 chromosomes as well as on 37 unplaced scaffolds. 203 of the 224 different DMRs (91%) overlapped genes or 5kb flanking regions, among which 52% were directly located in gene bodies, promoter or TSS sequences (35.3% in gene bodies and 47.3% in promoters/TSS) suggesting their putative functional role in gene expression and regulation and/or splicing events.
Following the procedure previously described, GO analyses on the pooled genes list revealed an overrepresentation of modules associated with the nervous system (GO:2000300, regulation of synaptic vesicle exocytosis; GO:0050804 modulation of synaptic transmission), immunity (GO:005728, negative regulation of inflammatory response; GO:0050852, T cell receptor signaling pathway), metabolic activity (GO:006816 calcium ion transport; GO:0055072, iron ion homeostasis, GO:0043087, regulation of GTPase activity), behaviour (GO:0007626, locomotory behavior) and endocrine processes (GO:0044060: regulation of endocrine process). All enriched GO are presented in Table 2 and SI Figure S5.
We also searched for DMRs between sexes, following the same procedure. We identified 206 DMRs associated with sex, of which 58 for Barcelona, 81 for Montpellier and 99 for Warsaw. Warsaw presented significantly more hyper than hypomethylated DMRs (χ2=5.878, P= 0.015), but it was not the case for Barcelona (χ2=0, P=1) nor Montpellier (χ2=1.25, P=0.264). DMRs were distributed on 29 chromosomes and 35 unplaced scaffolds. On a total of 206 DMRs, 181 (57.3%) were on genes or in a 5kb upstream/downstream region around genes. GO analyses revealed enrichment of genes involved in development, growth and morphogenesis, among others (see detailed enriched GO Table 3, Figure 6 & S6).
Almost twice more sex DMRs were shared between locations (11,7%) than between habitats (6.25%, see Figure S4 A & B ; z-test: χ²-squared  =  3.885, P  =  0.049). When taking into account the direction of methylation difference, 7 sex DMRs were shared between at least two cities (9.7%) which was three times more than for habitat DMRs (3.1%; z-test: X-squared  =  7.904, P  =  0.005).