# Disentangling group specific QTL allele effects from genetic background epistasis using admixed individuals in GWAS: an application to maize flowering

Simon Rio, Tristan Mary-Huard, Laurence Moreau, Cyril Bauland, Carine Palaffre, Delphine Madur, Valerie Combes, Alain Charcosset

► To cite this version:

Simon Rio, Tristan Mary-Huard, Laurence Moreau, Cyril Bauland, Carine Palaffre, et al. Disentangling group specific QTL allele effects from genetic background epistasis using admixed individuals in GWAS: an application to maize flowering. 2019. [hal-02789968](https://hal.inrae.fr/hal-02789968)

HAL Id: hal-02789968

https://hal.inrae.fr/hal-02789968

Preprint submitted on 5 Jun 2020

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Disentangling group specific QTL allele effects from genetic background epistasis using admixed individuals in GWAS: an application to maize flowering

Rio Simon¹, Mary-Huard Tristan¹², Moreau Laurence¹, Bauland Cyril¹, Palaffre Carine³, Madur Delphine¹, Combes Valérie¹, Charcosset Alain¹*,

1 GQE - Le Moulon, INRA, Univ. Paris-Sud, CNRS, AgroParisTech, Université Paris-Saclay, 91190, Gif-sur-Yvette, France

2 MIA, INRA, AgroParisTech, Université Paris-Saclay, 75005, Paris, France

3 UE 0394 SMH, INRA, 2297 Route de l'INRA, 40390, Saint-Martin-de-Hinx, France

* alain.charcosset@inra.fr

## Abstract

When handling a structured population in association mapping, group-specific allele effects may be observed at quantitative trait loci (QTLs) for several reasons: (i) a different linkage disequilibrium (LD) between SNPs and QTLs across groups, (ii) group-specific genetic mutations in QTL regions, and/or (iii) epistatic interactions between QTLs and other loci that have differentiated allele frequencies between groups. We present here a new genome-wide association (GWAS) approach to identify QTLs exhibiting such group-specific allele effects. We developed genetic materials including admixed progeny from different genetic groups with known genome-wide ancestries (local admixture). A dedicated statistical methodology was developed to analyze pure and admixed individuals jointly, allowing one to disentangle the factors causing the heterogeneity of allele effects across groups. This approach was applied to maize by developing an inbred "Flint-Dent" panel including admixed individuals that was evaluated for flowering time. Several associations were detected revealing a wide range of configurations of allele effects, both at known flowering QTLs (*Vgt1*, *Vgt2* and *Vgt3*) and new loci. We found several QTLs whose effect depended on the group ancestry of alleles while others interacted with the genetic background. The existence of directional epistasis was highlighted by comparing admixed with pure individuals and was consistent with epistatic interactions identified at the level of QTLs. Our GWAS approach provides useful information on the stability of QTL effects across genetic groups and can be applied to a wide range of species.

## Author summary

Identification of genomic regions involved in genetic architecture of traits has become commonplace in quantitative genetics studies. Genetic structure is a common feature in human, animal and plant species and most current methods target genomic regions whose effects on traits are conserved between genetic groups. However, a heterogeneity of allele effects may be observed due to different factors: a group-specific correlation between the alleles of the tagged marker and those of the causal variant, a group-specific mutation at the causal variant or an epistatic interaction between the causal variant and the genetic background. We propose a new method adapted to
structured populations including admixed individuals, which aims to identify these
genomic regions and to unravel the previous factors. method was applied to a maize
inbred diversity panel including lines from the dent and the flint genetic groups, as well
as admixed lines, evaluated for flowering time. Several genomic regions were detected
with various configurations of allele effects, with evidence of epistatic interactions
between some of the loci and the genetic background.

Introduction

Quantitative traits are genetically determined by numerous regions of the genome, also known as quantitative trait loci (QTLs). The advent of high density genotyping of single nucleotide polymorphism (SNPs) has opened the way to the identification of QTLs in diversity panels. These studies, referred to as genome-wide association studies (GWAS), use the linkage disequilibrium (LD) between the SNPs and the QTLs underlying the traits of interest. The panels evaluated in GWAS often include sets of individuals with complex pedigrees or genetic structure [1]. The latter is a common feature in human, animal and plant species and arises when groups of individuals cease to mate with each other and start to be subjected to different evolutionary forces.

Applying GWAS in a diversity panel including individuals from different groups raises the issue of spurious associations. The stratification of a population into genetic groups generates LD between loci that are differentiated between groups but not necessarily genetically linked. When a given trait is characterized by contrasted group-specific means, all these SNPs will correlate to it and may be detected as false positives. An efficient control of these spurious associations can be done by taking structure and kinship into account in the statistical model [1,2]. This procedure will however limit the statistical power at differentiated SNPs, making them difficult to detect in multi-group GWAS, especially in case of rare alleles [3].

In a structured population, group-specific allele effects can be observed at SNPs, and
testing an overall effect using a standard GWAS model may not be effective if the QTL
effect is of opposite sign in the different groups. Such effects can result from group
differences in LD between SNPs and QTLs across genetic groups. A different LD extent
or linkage phase between linked loci can be explained by specific dynamics of population
size such as bottlenecks or expansions [4,5]. Such patterns of LD were identified in
numerous species including human [6,7], dairy and beef cattle [8,9], pig [10], wheat [11]
and maize [12–15]. A genetic mutation appearing in a QTL region may also lead to
group-specific allele effects if it occurred in a founder specific of the genetic group.
Several Mendelian syndromes of obesity where shown to result from mutation within
specific ethnicities in human [16]. Another possibility consists in QTLs interacting with
other loci that have differentiated allele frequencies between groups (i.e. interacting
with the genetic background). In human, this possibility was discussed for a candidate
gene associated with a higher risk of myocardial infarction in African American than in
European populations [17,18]. Another example is a SNP in the promoter region of
*HNF4A* gene which was associated with a higher risk of developing type 2 diabetes in
Askenazi compared to United Kingdom populations [19]. This locus was later proven to
be interacting with another gene in the Askenazi population [20]. In maize, evidences of
QTLs with group-specific allele effects can also be found, even though the cause of these
differences remains unclear. The presence of allelic series has been demonstrated for
QTLs associated with flowering time, including *Vgt1* [21]. A QTL with group-specific
allele effects was also identified in a maize diversity panel for a phenology trait [22].
More generally, studying the stability of QTL allele effects across genetic backgrounds is
an important issue. In human, it determines the ability of a genetic marker to predict
the predisposition of an individual to develop a genetic disease across ethnic groups. In
plant or animal breeding, it conditions the success of introgressing a favorable allele coming from a source of diversity into an elite genetic material.

Different GWAS strategies were adopted to address this issue depending on the species. In human, GWAS mostly focused on a specific genetic group, and these group-specific studies were compared later through meta-analyses [23,24]. Some of these meta-analyses revealed highly conserved effects between populations [25, 26] while other put in evidence more differences [27]. In dairy cattle, the first GWAS studies focused on a specific breed [28–30]. More recently, multi-breed GWAS were conducted to refine QTLs locations by taking advantage of the low LD extent observed in such composite populations [31–33]. In maize, the possibility to use seeds from different origins and generations led geneticists to assemble GWAS panel with a broad range of genetic materials [34–36]. These panels often include a limited proportion of admixed individuals that were derived from crosses between individuals from different genetic groups. The genomes of these admixed individuals consist in mosaics of fragments with different ancestries. Admixture events are a common feature in living species and can contribute to the successful colonization of new environments [37, 38]. In plants, innovative admixed genetic materials were created to enable high statistical power of QTL detection along with a wide spectrum of genetic diversity studied, such as nested association mapping (NAM) [39] or multi-parent advanced generation inter-cross (MAGIC) [40]. Both NAM and MAGIC populations are of great interest to study the stability of QTL effects in a wide range of genetic backgrounds. However, they generally include a limited number of founders and do not address the stability of QTL allele effects across genetic groups.

This study aimed at evaluating the interest of producing admixed individuals, derived from a large set of parents, in order to decipher the genetic architecture of a trait using innovative GWAS models. The objectives were (i) to demonstrate the interest of multi-group analyses to identify new QTLs, (ii) to highlight the interest of applying multi-group GWAS models to identify group-specific allele effects at QTLs and (iii) to show how admixed individuals can help to disentangle the factors causing the heterogeneity of allele effects across groups: local genomic differences or epistatic interactions between QTLs and the genetic background. To our knowledge, no method has been proposed in the literature to address the last objective. This method was applied to a maize inbred population evaluated for flowering traits, including dent, flint and admixed lines. Maize flowering time is an interesting trait to analyze in quantitative genetics studies. It is considered as a major adaptive trait by tailoring vegetative and reproductive growth phases to local environmental conditions. Admixed individuals were also used to investigate the existence of directional epistasis using a test based on the mean of admixed individuals relative to that of their progenitors.

# Materials and methods

## Genetic material and genotypic data

Genetic material consisted in a panel of 970 maize inbred lines assembled within the "Amaizing" project. It gathered 300 dent lines, 304 flint lines and 366 admixed doubled haploids, further referred to as admixed lines. The dent lines were those included in the "Amaizing Dent" panel [41] and the flint lines were those included in the "CF-Flint" panel [15]. The dent and flint lines aimed at representing the diversity of their respective heterotic group used in European breeding and included several breeding generations. The admixed lines were derived from 206 hybrids between flint and dent lines, mated according to a sparse factorial design (Fig 1), followed by in situ gynogenesis [42] to produce fixed admixed inbred lines. Each dent or flint line was involved in 0 to 11
hybrids (1.21 in average), each leading to 1 to 4 admixed lines (1.77 in average). In total, 171 dent lines and 172 flint lines were involved as parents of admixed lines.

94

95

<figure><img src="image_2.png" /><figcaption>**Fig 1.** Diagram of admixed lines production from hybrids obtained by mating dent and flint lines according to a sparse factorial design</figcaption></figure>

All the flint and dent lines were genotyped using the 600K Affymetrix Maize Genotyping Array [43]. Residual heterozygous data was treated as missing and all missing values were imputed independently within each group using Beagle v.3.3.2 and default parameters (Browning and Browning 2009). The admixed lines were genotyped with a 15K chip provided by the private company Limagrain which included a reduced set of SNPs from the 50K Illumina MaizeSNP50 BeadChip [44]. Eight check lines were included in both datasets to standardize the allele coding (0/1) on the common SNPs (around 9,000). The following procedure to impute admixed genotypes up to 600K SNPs is illustrated in S1 Fig. The positions of recombination breakpoints and the parental origin of the alleles for admixed lines were determined with these common SNPs. A smoothing of parental allele origins was performed for the few SNPs indicating discordant information with respect to the chromosome block in which they were located. In this case, we considered the underlying genotypic datapoint as missing. Parental origins of alleles in admixed lines were imputed up to 600K using adjacent SNP information. If a set of SNPs to be imputed was located within a recombination interval, the new position of the breakpoint was positioned at half of that ordered set, according to the physical position of the SNPs along the chromosome. Alleles at SNPs were then imputed based on their origin using parental genotypic data. The MITE associated with the flowering QTL *Vgt1* [45, 46] was also genotyped for all the individuals (0: absence, 1: presence). There was a total of 482,013 polymorphic SNPs in this dataset, for which we had information for each individual concerning the SNP allele (0/1), its ancestry (dent/flint) and the genetic background (dent/flint/admixed) in which it was observed.

The dent genome proportion of the admixed lines ranged from 0.16 to 0.86 with a mean equal to 0.51 (S2 Fig). Possible selection biases were studied along the genome by comparing the observed allele frequencies with the expected allele frequencies given the pedigree. No major pattern was observed, suggesting no or minor selection biases among the admixed lines (S3 Fig). A PCoA was performed on genetic distances computed as $D_{i,j} = 1 - K_{i,j}$, with $K_{i,j}$ being the kinship coefficient between lines *i* and *j* computed following Eq (2, see below). The flint and dent lines are clearly distinguished on the two

96

97

98

99

100

101

102

103

104

105

106

107

108

109

110

111

112

113

114

115

116

117

118

119

120

121

122

123

124
principal components, with a small overlapping region in the center of the graph, while
the admixed lines fill the genetic space between the two groups (Fig 2).

125

126

<figure><img src="image_3.png" /><figcaption>**Fig 2.** PCoA on genetic distances with coloration of individuals depending on their type: dent, flint or admixed lines</figcaption></figure>

LD between pairs of loci was estimated separately in the dent and the flint datasets using two estimators. The first was the standard measure of LD, computed as the square correlation between pairs of loci $r^2$. The second was the estimator $r_K^2$ proposed by [47], accounting for relatedness estimated using Eq (2). We only considered SNPs for which at least ten individuals carried the minor allele in both dent and flint datasets. LD extent was compared between groups using both estimators. A sliding window of 1Kbp up to 2 Mbp was used to group pairs of loci with similar physical distances. The average LD was computed within each group-specific windows and revealed a higher LD extent in the dent than in the flint genetic group (S4 Fig), which was consistent with previous studies [12–15]. As suggested by [8], the persistence of LD linkage phases across flint and dent genetic groups was evaluated by computing the correlation between the $r$ (and $r_K$) estimated in each group using a sliding window of 1Kbp up to 2 Mbp. We also studied the consistency of marker phases between groups by computing, for each LD estimator, the correlation between their signs in the two groups. LD phases were very consistent over short physical distances but began to diverge dramatically when the loci were distant by more than 100-200 Kbp (S5 Fig).

127

128

129

130

131

132

133

134

135

136

137

138

139

140

141

142

Phenotypic data

143

All the lines were evaluated *per se* at Saint-Martin-de-Hinx (France) in 2015 and 2016 for male flowering (MF) and female flowering (FF), in calendar days after sowing. Each plot consisted in a row of 25 plants. MF and FF were measured as a median value within the whole plot. Each trial was a latinized alpha design where every line was

144

145

146

147
evaluated two times on average. Field trials were divided into blocks of 36 plots each. To avoid competition between genetic backgrounds, dent, flint and admixed lines were sown in different blocks. Three check individuals were repeated in all blocks (B73, F353 and UH007).

Variance components were estimated using model:

$$ 
\begin{align*}
Y_{jklrc} &= \mu + \beta_j + \alpha_k + G_{lk} + (G \times \beta)_{lkj} + X_{rj} + Z_{cj} + E_{jklrc} \\
G_{lk} &\sim \mathcal{N}(0, \sigma_{G_k}^2) \text{ independent} \\
(G \times \beta)_{lkj} &\sim \mathcal{N}(0, \sigma_{(G \times \beta)_{jk}}^2) \text{ independent} \\
E_{jklrc} &\sim \mathcal{N}(0, \sigma_{E_j}^2) \text{ independent,}
\end{align*}
$$

where $Y_{jklrc}$ is the phenotype, $\mu$ is the intercept, $\beta_j$ is the fixed effect of trial j, $\alpha_k$ is the fixed effect of genetic background k (dent, flint, admixed, or the different checks: B73, F353 and UH007), $G_{lk}$ is the random genotype effect of line l in genetic background k (not for checks) with $\sigma_{G_k}^2$ being the genotypic variance in genetic background k, $(G \times \beta)_{lkj}$ is the random Genotype x Environment (GxE) interaction of line l in genetic background k for trial j, with $\sigma_{(G \times \beta)_{jk}}^2$ being the GxE variance in the genetic background k for trial j, $E_{jklrc}$ is the error with $\sigma_{E_j}^2$ being the error variance for trial j, $X_{rj}$ and $Z_{cj}$ are the row and column random effects in trial j respectively, as defined by the field design. All random effects are independent of each other. The row and column effects were modeled as independent or using an autoregressive model (AR1), as determined based on the AIC criterion (S1 Table). Least squares means, further referred to as phenotypes, were computed over the whole design using the same model, with genotypes as fixed effects. Model parameters were estimated using ASReml-R and restricted maximum likelihood (ReML) [48].

## Global assessment of directional epistasis

This panel allowed us to test for the existence of directional epistasis, which refers to epistatic interactions that are biased toward high or low genetic values [49]. In the presence of directional epistatic interactions and provided no selection, we can expect the genetic mean of the admixed lines to be different from its expected value, obtained by considering only additive effects (S1 Appendix). The existence of directional epistasis was investigated using a test based on the comparison between the means of the progeny and the parental populations. The following model was applied on the joint dent, flint and admixed dataset:

$$ Y_{kl} = \mu + \alpha_k + G_{kl} + E_{kl} \quad (1) $$

where $Y_{kl}$ is the phenotype (least squares means) of the line among the N individuals of the sample, $\mu$ is the intercept, $\alpha_k$ is the genetic background effect with $k \in \{D, F, A\}$ for dent, flint and admixed genetic background respectively. $G_{kl}$ is the random genetic value of the line where $g$ is the vector of genetic values with $g \sim \mathcal{N}(0, K\sigma_G^2)$, $K$ is the kinship matrix computed following Eq (2) using allele frequencies estimated on the joint dent, flint and admixed dataset, $\sigma_G^2$ is the genetic variance, $E_{kl}$ is the residual error of the line where $e$ is the vector of residuals with $e \sim \mathcal{N}(0, I\sigma_E^2)$, $I$ is the identity matrix and $\sigma_E^2$ is the residual variance. For each trait, the linear combination $H_0: \frac{1}{2}(\alpha_D + \alpha_F) - \alpha_A = 0$ was tested to identify directional epistasis.

The kinship between individuals i and j, $K_{ij}$, was computed following [50]:

$$ K_{ij} = \frac{\sum_{m=1}^{M} (W_{im} - f_m)(W_{jm} - f_m)}{\sum_{m=1}^{M} f_m(1 - f_m)} \quad (2) $$
where $W_{im}$ is the genotype of individual $i$ at locus $m$ coded 0/1 and $f_m$ is the frequency of allele 1 at locus $m$.

## GWAS models

In this study, three GWAS models were applied to different population samples (Table 1). The GWAS strategies were (i) to analyze dent and flint lines separately using a standard GWAS model $M_1$, (ii) to analyze dent and flint lines jointly using a GWAS model $M_2$ accounting for allele ancestry (confounded with the genetic background) and (iii) to analyze dent, flint and admixed lines in a GWAS model $M_3$ accounting for both allele ancestry and the genetic background of the individuals. All models aimed at detecting a SNP effect, defined as a contrast effect between alleles 0 and 1 at a given SNP.

<table><caption>**Table 1.** Population sample to which each GWAS model was applied with the corresponding number of SNPs conserved for the analysis (at least 10 individuals carrying the minor allelic state)</caption><thead><tr><th></th><th>Dent</th><th>Flint</th><th>Dent + Flint</th><th>Dent + Flint + Admixed</th></tr></thead><tbody><tr><th>M<sub>1</sub></th><td>✓ (247,759)</td><td>✓ (282,278)</td><td>✗</td><td>✗</td></tr><tr><th>M<sub>2</sub></th><td>-</td><td>-</td><td>✓ (288,093)</td><td>✗</td></tr><tr><th>M<sub>3</sub></th><td>-</td><td>-</td><td>-</td><td>✓ (256,951)</td></tr></tbody></table>

✓: model was applied to the sample

✗: model was not applied to the sample but can theoretically be, provided the addition of a genetic background effect

-: model cannot be applied to the sample or would simplify into another model

Note that the number of SNP in multi-group GWAS ($M_2$, $M_3$) is higher than the minimum number of SNPs in single group GWAS ($M_1$ (Dent)). SNPs carrying redundant information within a single group were indeed reduced to a single SNP for M1 and may no longer carry redundant information when datasets are pooled ($M_2$, $M_3$)

### Standard GWAS model $M_1$

The first GWAS model $M_1$ [1] was applied separately to the dent and flint datasets. For each SNP among the $M$ loci, one has:

$$Y_{il} = \mu + \beta_{i}^{m} + G_{il} + E_{il}$$

where $\beta_i^m$ is the effect of the SNP allele $i$ at locus $m$ (Table 2). All other terms are identical to those described Eq (1), and the kinship was computed following Eq (2) using allele frequencies estimated for each dataset. The existence of a SNP effect was tested using hypothesis $H_0: \Delta^m = \beta_1^m - \beta_0^m = 0$.

### Multi-group GWAS model $M_2$

We applied a multi-group GWAS model $M_2$ jointly to the flint and dent datasets, specifying the allele ancestry (confounded with the genetic background). For a given SNP $m$, one has:

$$Y_{ijl} = \mu + \beta_{ij}^{m} + G_{ijl} + E_{ijl}$$

where $\beta_{ij}^m$ is the effect of the SNP allele $i$ with ancestry $j$ at locus $m$, as defined in Table 2. All other terms are identical to those described Eq. (1), and the kinship was computed following Eq (2) using allele frequencies estimated on the joint dent and flint dataset. At a given SNP, the following hypotheses were tested:
<table><caption>**Table 2.** Allelic states observed in each GWAS model, resulting from a combination of SNP alleles, their ancestry and the genetic background in which they are observed.</caption><thead><tr><th></th><th>SNP</th><th>Ancestry</th><th>Genetic background</th><th>Allelic states</th></tr></thead><tbody><tr><td>M<sub>1</sub></td><td>{0, 1}</td><td>-</td><td>-</td><td>{0, 1}</td></tr><tr><td>M<sub>2</sub></td><td>{0, 1}</td><td>{D, F}<sup>a</sup></td><td>-</td><td>{0D, 1D, 0F, 1F}</td></tr><tr><td>M<sub>3</sub></td><td>{0, 1}</td><td>{D, F}</td><td>{D, A, F}</td><td>{0DD, 1DD, 0DA, 1DA, 0FA, 1FA, 0FF, 1FF}</td></tr></tbody></table>

0 : SNP reference allele

1 : SNP alternative allele

D : Dent ancestry or genetic background

F : Flint ancestry or genetic background

A : Admixed genetic background

a confounded with the genetic background

$$ \bullet H_0 : \Delta_D^m = \beta_{1D}^m - \beta_{0D}^m = 0 $$

$$ \bullet H_0 : \Delta_F^m = \beta_{1F}^m - \beta_{0F}^m = 0 $$

$$ \bullet H_0 : \Delta_{D+F}^m = \Delta_D^m + \Delta_F^m = 0 $$

$$ \bullet H_0 : \Delta_{D-F}^m = \Delta_D^m - \Delta_F^m = 0 $$

Hypotheses $\Delta_D^m$ and $\Delta_F^m$ test the existence of a dent and a flint SNP effect respectively. Hypothesis $\Delta_{D+F}^m$ tests for a general SNP effect while $\Delta_{D-F}^m$ tests for a divergent SNP effect between the dent and flint ancestries

**Multi-group GWAS model M<sub>3</sub>**

We applied a multi-group GWAS model M<sub>3</sub> jointly to the flint, dent and admixed datasets, specifying the allele ancestry and the genetic background of the individual. For a given SNP *m*, one has:

$$ Y_{ijkl} = \mu + \beta_{ijk}^{m} + G_{ijkl} + E_{ijkl} $$

where $\beta_{ijk}^m$ is the effect of the SNP allele *i* with ancestry *j* at locus *m* in genetic background *k*, as defined in Table 2. All other terms are identical to those described in Eq (1), and the kinship was computed following Eq (2) using allele frequencies estimated on the joint dent, flint and admixed dataset. At a given SNP, 16 hypotheses were tested (Table 3). Hypotheses referred to as "simple" ($\Delta_{DD}^m$, $\Delta_{DA}^m$, $\Delta_{FA}^m$ and $\Delta_{FF}^m$) were tested to identify QTLs with a significant SNP effect for each combination of ancestries and genetic backgrounds. For instance, $\Delta_{DD}^m$ tests whether a dent SNP effect (differential effect between alleles 0 and 1 of dent ancestry) exists in the dent genetic background. Hypotheses referred to as "general" ($\Delta_{FA+FF}^m$, $\Delta_{DD+DA}^m$, $\Delta_{DA+FA}^m$, $\Delta_{DD+FF}^m$ and, $\Delta_{DD+DA+FA+FF}^m$) were used to identify QTLs with a mean SNP effect over ancestries and genetic backgrounds. For instance, $\Delta_{FA+FF}^m$ tests for a general flint SNP effect in the flint and the admixed genetic backgrounds and $\Delta_{DD+DA+FA+FF}^m$ tests for a general SNP effect over ancestries and genetic backgrounds. Hypotheses referred to as "divergent" ($\Delta_{DA-FA}^m$, $\Delta_{DD-DA}^m$, $\Delta_{FA-FF}^m$, $\Delta_{DD-FF}^m$, $\Delta_{DA-FF}^m$, $\Delta_{DD-FA}^m$, $\Delta_{(DD+DA)-(FA+FF)}^m$, $\Delta_{(DD+FF)-(DA+FA)}^m$, $\Delta_{(DD-DA)-(FF-FA)}^m$) were tested to identify QTLs with a contrasted SNP effect between ancestries and/or genetic backgrounds. For instance, $\Delta_{DD-DA}^m$ tests for a divergent dent SNP effect between the dent and the admixed genetic backgrounds, which amounts to testing an epistatic interaction between the SNP and the genetic background (see S2 Appendix for details).

On a biological standpoint, a QTL with contrasted SNP effects between groups can be caused by (i) a local genomic difference due to a group-specific genetic mutation
<table><caption>**Table 3.** Linear combination tested with $M_3$ compared to hypotheses tested using other GWAS models ($M_1$ and $M_2$).</caption><thead><tr><th></th><th>Type</th><th>Δ<sup>m</sup><sub>DD</sub><sup>a</sup></th><th>Δ<sup>m</sup><sub>DA</sub><sup>b</sup></th><th>Δ<sup>m</sup><sub>FA</sub><sup>c</sup></th><th>Δ<sup>m</sup><sub>FF</sub><sup>d</sup></th><th>M<sub>1</sub></th><th>M<sub>2</sub></th></tr></thead><tbody><tr><td>Δ<sup>m</sup><sub>DD</sub></td><td>simple</td><td>+1</td><td>0</td><td>0</td><td>0</td><td>✓</td><td>✓</td></tr><tr><td>Δ<sup>m</sup><sub>DA</sub></td><td>simple</td><td>0</td><td>+1</td><td>0</td><td>0</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>FA</sub></td><td>simple</td><td>0</td><td>0</td><td>+1</td><td>0</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>FF</sub></td><td>simple</td><td>0</td><td>0</td><td>0</td><td>+1</td><td>✓</td><td>✓</td></tr><tr><td>Δ<sup>m</sup><sub>FA+FF</sub></td><td>general</td><td>0</td><td>0</td><td>+1</td><td>+1</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>DD+DA</sub></td><td>general</td><td>+1</td><td>+1</td><td>0</td><td>0</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>DA+FA</sub></td><td>general</td><td>0</td><td>+1</td><td>+1</td><td>0</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>DD+FF</sub></td><td>general</td><td>+1</td><td>0</td><td>0</td><td>+1</td><td>-</td><td>✓</td></tr><tr><td>Δ<sup>m</sup><sub>DD+DA+FA+FF</sub></td><td>general</td><td>+1</td><td>+1</td><td>+1</td><td>+1</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>DA-FA</sub></td><td>divergent</td><td>0</td><td>+1</td><td>-1</td><td>0</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>DD-DA</sub></td><td>divergent</td><td>+1</td><td>-1</td><td>0</td><td>0</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>FA-FF</sub></td><td>divergent</td><td>0</td><td>0</td><td>+1</td><td>-1</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>DD-FF</sub></td><td>divergent</td><td>+1</td><td>0</td><td>0</td><td>-1</td><td>-</td><td>✓</td></tr><tr><td>Δ<sup>m</sup><sub>(DD+DA)-(FA+FF)</sub></td><td>divergent</td><td>+1</td><td>+1</td><td>-1</td><td>-1</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>(DD+FF)-(DA+FA)</sub></td><td>divergent</td><td>+1</td><td>-1</td><td>-1</td><td>+1</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>(DD-DA)-(FF-FA)</sub></td><td>divergent</td><td>+1</td><td>-1</td><td>+1</td><td>-1</td><td>-</td><td>-</td></tr></tbody></table>

✓: hypothesis also tested using the corresponding GWAS model

-: hypothesis not tested using the corresponding GWAS model

and/or to group differences in LD or (ii) an interaction with the genetic background. Under the first hypothesis, one expects that the effect of a SNP depends on its ancestry but not on the genetic background (admixed or pure, see Fig 3-a). Under the second hypothesis, we expect a SNP effect, for a given ancestry, to vary depending on the genetic background. One example would be a QTL with a strong SNP effect in a dent genetic background, but none in the flint genetic background, while the SNP effects would be of intermediate size for alleles of both ancestries in the admixed genetic background (see Fig 3-b). Note that other complex configurations are possible, justifying the inclusion of all tests in the analysis.

For the three GWAS models, a SNP was discarded if its minor allelic state (Table 2) was carried by less than 10 individuals, or if it carried a redundant information with another SNP. Model parameters were estimated using ReML and the linear combinations of fixed effects were tested using Wald tests, both implemented in the R-package MM4LMM [51]. The false discovery rate (FDR) was controlled by applying the procedure of [52] jointly to the whole set of tests defined by each GWAS strategy, and repeatedly for each trait. For a given hypothesis tested, significant SNPs were clustered into QTLs if they were located within a physical window of 3 Mbp, leading to a LD below 0.05 between markers of different QTLs.

<footnote>a $Δ_{DD}^m = β_{1DD}^m - β_{0DD}^m$</footnote>

<footnote>b $Δ_{DA}^m = β_{1DA}^m - β_{0DA}^m$</footnote>

<footnote>c $Δ_{FA}^m = β_{1FA}^m - β_{0FA}^m$</footnote>

<footnote>d $Δ_{FF}^m = β_{1FF}^m - β_{0FF}^m$</footnote>
<figure><img src="image_4.png" /><figcaption>**Fig 3.** Schematic of allele effects when divergent SNP effects are observed between groups, depending on the biological hypothesis: (a) local genomic difference between groups and (b) allele effects interacting with the genetic background. The denomination of the allelic states on the x-axis include the SNP allele (0/1), its ancestry (D/F) and the genetic background in which it is observed (D/A/F), as presented in Table 2.</figcaption></figure>

# Results

## Phenotypic analysis and directional epistasis

We observed a substantial phenotypic variability within the dent, flint and admixed genetic backgrounds. The variance components estimated in the phenotypic analysis were summarized in S1 Table. Similar trends were observed for both MF and FF. The admixed genotypic variance was lower than the dent and flint genotypic variances, which were themselves comparable. GxE variances were limited and the broad sense heritabilities were high for each genetic background, ranging from 0.88 in the admixed lines to 0.96 in the dent and flint lines for both MF and FF.

The presence of admixed lines allowed us to test the existence of directional epistasis which was significant for both MF and FF (Table. 4). The mean of admixed lines estimated using a model accounting for relatedness differed significantly from the one expected without directional epistatic interactions. On average, admixed lines flowered as late as dent lines while the flint lines flowered earlier.

<table><caption>**Table 4.** Test for directional epistasis with group-specific means estimated by the model (Eq.1) and the p-value (pval) of the directional epistatic deviation</caption><thead><tr><th></th><th>Dent</th><th>Flint</th><th>Admixed</th><th>pval</th></tr></thead><tbody><tr><td>MF</td><td>68.26</td><td>66.26</td><td>68.44</td><td>3.14 10<sup>-10</sup> ***</td></tr><tr><td>FF</td><td>69.84</td><td>67.87</td><td>70.16</td><td>2.05 10<sup>-11</sup> ***</td></tr></tbody></table>

## Associations detected and comparison of GWAS strategies

For each GWAS model, two levels of FDR were used: 5% and 20% to declare a SNP as significantly associated. The number of significant SNPs detected and the corresponding number of QTLs were summarized in Table 5 for both traits. The location of QTLs detected using a FDR of 20% was represented along the genome in Fig

<footnote>*** : pval < 10⁻³ ; ** : < 10⁻³ pval < 10⁻² ; * : < 10⁻² pval < 5 × 10⁻²</footnote>
4 for MF and in S6 Fig for FF. All associations were listed in S2 Table and S3 Table. Note that major QTLs detected by a model (e.g. M₁) may be discarded with another model (e.g. M₃) because of the filtering on allele frequencies.

272

273

274

<table><caption>**Table 5.** Number of SNPs associated with each trait, depending on the GWAS strategy, using a FDR of 5% and 20%. The number of corresponding QTLs is also indicated</caption><thead><tr><th rowspan="3"></th><th colspan="4">MF</th><th colspan="4">FF</th></tr><tr><th colspan="2">5%</th><th colspan="2">20%</th><th colspan="2">5%</th><th colspan="2">20%</th></tr><tr><th>SNP</th><th>QTL</th><th>SNP</th><th>QTL</th><th>SNP</th><th>QTL</th><th>SNP</th><th>QTL</th></tr></thead><tbody><tr><td>M<sub>1</sub><sup>a</sup></td><td>7</td><td>2</td><td>56</td><td>24</td><td>8</td><td>3</td><td>38</td><td>14</td></tr><tr><td>Δ<sup>m</sup> (Dent)</td><td>4</td><td>1</td><td>35</td><td>12</td><td>4</td><td>1</td><td>22</td><td>6</td></tr><tr><td>Δ<sup>m</sup> (Flint)</td><td>3</td><td>1</td><td>21</td><td>13</td><td>4</td><td>2</td><td>16</td><td>8</td></tr><tr><td>M<sub>2</sub><sup>a</sup></td><td>4</td><td>1</td><td>57</td><td>22</td><td>4</td><td>1</td><td>7</td><td>3</td></tr><tr><td>Δ<sup>m</sup><sub>D</sub></td><td>4</td><td>1</td><td>37</td><td>9</td><td>4</td><td>1</td><td>4</td><td>1</td></tr><tr><td>Δ<sup>m</sup><sub>F</sub></td><td>-</td><td>-</td><td>3</td><td>2</td><td>-</td><td>-</td><td>2</td><td>1</td></tr><tr><td>Δ<sup>m</sup><sub>D+F</sub></td><td>1</td><td>1</td><td>11</td><td>7</td><td>2</td><td>1</td><td>2</td><td>1</td></tr><tr><td>Δ<sup>m</sup><sub>D-F</sub></td><td>-</td><td>-</td><td>18</td><td>9</td><td>-</td><td>-</td><td>1</td><td>1</td></tr><tr><td>M<sub>3</sub><sup>a</sup></td><td>7</td><td>3</td><td>116</td><td>41</td><td>6</td><td>2</td><td>11</td><td>6</td></tr><tr><td>Δ<sup>m</sup><sub>DD</sub></td><td>4</td><td>1</td><td>32</td><td>8</td><td>4</td><td>1</td><td>4</td><td>1</td></tr><tr><td>Δ<sup>m</sup><sub>DA</sub></td><td>-</td><td>-</td><td>1</td><td>1</td><td>-</td><td>-</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>FA</sub></td><td>2</td><td>1</td><td>10</td><td>2</td><td>-</td><td>-</td><td>1</td><td>1</td></tr><tr><td>Δ<sup>m</sup><sub>FF</sub></td><td>-</td><td>-</td><td>1</td><td>1</td><td>-</td><td>-</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>FA+FF</sub></td><td>-</td><td>-</td><td>4</td><td>4</td><td>-</td><td>-</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>DA+DD</sub></td><td>-</td><td>-</td><td>10</td><td>4</td><td>-</td><td>-</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>DA+FA</sub></td><td>-</td><td>-</td><td>11</td><td>5</td><td>-</td><td>-</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>DD+FF</sub></td><td>2</td><td>2</td><td>34</td><td>12</td><td>2</td><td>2</td><td>6</td><td>2</td></tr><tr><td>Δ<sup>m</sup><sub>DD+DA+FA+FF</sub></td><td>-</td><td>-</td><td>19</td><td>6</td><td>1</td><td>1</td><td>2</td><td>1</td></tr><tr><td>Δ<sup>m</sup><sub>DA-FA</sub></td><td>-</td><td>-</td><td>2</td><td>2</td><td>-</td><td>-</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>DD-DA<sup>b</sup></sub></td><td>-</td><td>-</td><td>4</td><td>4</td><td>-</td><td>-</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>FA-FF</sub></td><td>-</td><td>-</td><td>1</td><td>1</td><td>-</td><td>-</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>DD-FF</sub></td><td>-</td><td>-</td><td>15</td><td>5</td><td>-</td><td>-</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>(DD+DA)-(FA+FF)</sub></td><td>-</td><td>-</td><td>5</td><td>3</td><td>-</td><td>-</td><td>-</td><td>-</td></tr><tr><td>Δ<sup>m</sup><sub>(DD+FF)-(DA+FA)</sub><sup>b</sup></td><td>-</td><td>-</td><td>5</td><td>4</td><td>-</td><td>-</td><td>2</td><td>2</td></tr><tr><td>Δ<sup>m</sup><sub>(DD-DA)-(FF-FA)</sub></td><td>-</td><td>-</td><td>2</td><td>2</td><td>-</td><td>-</td><td>-</td><td>-</td></tr></tbody></table>

a number of SNPs detected over the set of tests (a given SNP can be detected using different tests)

b hypothesis testing an interaction between the QTL and the genetic background

First, a standard GWAS model M₁ was applied separately to the dent and the flint datasets. Based on a 20% FDR, 35 SNPs were associated with MF in the dent dataset while 21 SNPs were associated in the flint dataset. These SNPs can be clustered into 12 QTLs in the dent dataset and into 13 QTLs in the flint dataset. Interestingly, none of these SNPs were detected in both datasets and they only pointed to one common QTL between datasets, which was located in the vicinity of Vgt2 on chromosome 8 [14].

Secondly, dent and flint datasets were analyzed jointly using model M₂, which takes into account the dent or flint ancestry of the allele. Note that the allele ancestry is confounded with the genetic background in this model. Based on a 20% FDR, 57 SNPs were associated with MF and were significant for Δ<sup>m</sup><sub>D</sub> (37 SNPs), Δ<sup>m</sup><sub>F</sub> (3 SNPs), Δ<sup>m</sup><sub>D+F</sub> (11 SNPs) and Δ<sup>m</sup><sub>D-F</sub> (18 SNPs). Note that some SNPs displayed more than one significant test, which explains why the total number of SNPs over the four tests did not sum to 57. These SNPs can be clustered into 22 QTLs that were significant for Δ<sup>m</sup><sub>D</sub> (9

275

276

277

278

279

280

281

282

283

284

285

286

287
<figure><img src="image_5.png" /><figcaption>**Fig 4.** Position of QTLs detected by each GWAS strategy for MF using a FDR of 20%. The size of the grey dots is proportional to the -log₁₀(pval) of the test at the most significant SNP of the region. Red vertical lines and names below correspond to QTL discussed in section "Highlighted QTLs". Note that major QTLs detected by a model may be discarded with another model because of filtering on allele frequencies</figcaption></figure>

QTLs), $\Delta_F^m$ (2 QTLs), $\Delta_{D+F}^m$ (7 QTLs), $\Delta_{D-F}^m$ (9 QTLs). Note that some QTLs were already detected using $M_1$ such as the QTL located in the vicinity of *Vgt3* on chromosome 3 [53, 54] detected in the dent dataset. Other QTLs were specific to $M_1$ such as the QTL located on chromosome 2 detected in the flint dataset, or specific to $M_2$ such as the QTL located chromosome 5 detected using $\Delta_{D-F}^m$. Based on a 20% FDR, a similar number of QTLs was detected between $M_1$ and $M_2$ for MF, while more QTLs were detected using $M_1$ than $M_2$ for FF.

Finally, the dent, flint and admixed lines were analyzed jointly using model $M_3$ which distinguished the allele ancestry and the genetic background. The existence of a dent SNP effect was tested in the dent ($\Delta_{DD}^m$) and in the admixed genetic backgrounds ($\Delta_{DA}^m$), and similarly for the flint SNP effect ($\Delta_{FF}^m$ and $\Delta_{FA}^m$). Several hypotheses on general and divergent SNP effects were also tested between ancestries and genetic backgrounds (Table 3). Based on a 20% FDR, 116 SNPs were associated with MF and were significant for $\Delta_{DD}^m$ (32 SNPs), $\Delta_{DD+FF}^m$ (34 SNPs), $\Delta_{DD-DA}^m$ (4 SNPs),
$\Delta_{(DD+DA)-(FA+FF)}^m$ (5 SNPs) and others. These SNPs can be clustered into 41 QTLs that were significant for $\Delta_{DD}^m$ (8 QTLs), $\Delta_{DD+FF}^m$ (12 QTLs), $\Delta_{DD-DA}^m$ (4 QTLs), $\Delta_{(DD+DA)-(FA+FF)}^m$ (3 QTLs) and others. Note that some of the QTLs were already detected using $\mathbf{M}_1$ and $\mathbf{M}_2$ such as the QTL located in the vicinity of $Vgt3$ on chromosome 3, while several QTLs were specific to $\mathbf{M}_3$ such as the two QTLs detected in chromosome 2 using $\Delta_{FA}^m$. Several QTLs were detected as showing a divergent SNP effect, including hypotheses testing an interaction with the genetic background. Based on 5% and 20% FDRs, the number of QTLs detected with $\mathbf{M}_3$ was the highest for MF and intermediate between $\mathbf{M}_1$ and $\mathbf{M}_2$ for FF.

## Highlighted QTLs

Among the 41 QTLs detected for MF with $\mathbf{M}_3$, six QTLs were selected and studied in further details. The five first QTLs had (i) at least one significant test among $\mathbf{M}_3$ hypotheses based on a FDR of 20%, and (ii) a large frequency for each allele with a minimum of 30 lines carrying the minor allelic state ($QTL7.2$). Among them, one SNP was located in the vicinity of $Vgt2$ [14] and another in the vicinity of $Vgt3$ [53, 54]. In addition to these five QTLs, we also considered a MITE polymorphism known to be associated with $Vgt1$, a flowering QTL detected in several studies [21, 45, 46]. For all QTLs, information concerning their physical position along the genome, the frequency of each allelic state and their $-\log_{10}(\text{pval})$ at each test was summarized in Table 6. The distribution of the phenotypes is illustrated for each allele after adjusting the variation due to the polygenic background in Fig. 5, and their location along the genome is indicated by red vertical lines in Fig. 4. Other QTLs had interesting profiles, showing either group-specific allele effects conserved between ancestries or interactions with the genetic background, and are presented in S7 Fig and S4 Table.

The SNP matching $Vgt2$ region on chromosome 8 was detected as associated with MF (5% FDR) using $\Delta_{DD+FF}^m$ ($-\log_{10}(\text{pval})=7.15$) in $\mathbf{M}_3$. This QTL showed a conserved effect across ancestries and genetic backgrounds (Fig. 5-a). This observation was supported by a high $-\log_{10}(\text{pval})$ for tests related to a general SNP effect: $\Delta_{D+F}^m$ (5.25), $\Delta_{DD+DA}^m$ (5.35), $\Delta_{DA+FA}^m$ (3.20), $\Delta_{DD+FF}^m$ (7.15) and $\Delta_{DD+DA+FA+FF}^m$ (6.46), and a low $-\log_{10}(\text{pval})$ for tests related to divergent SNP effects (all below 1).

The SNP matching $Vgt3$ region on chromosome 3 was detected as associated with MF (5% FDR) using $\Delta_{DD}^m$ (10.53) in $\mathbf{M}_3$. This QTL showed a large effect in the dent genetic background, a medium effect in the admixed genetic background regardless of the allele ancestry and a small effect in the flint genetic background (Fig. 5-b). This observation was supported by a high $-\log_{10}(\text{pval})$ for the tests related to the dent SNP effect in the dent genetic background: $\Delta^m$ (Dent, 10.99), $\Delta_D^m$ (9.65) and $\Delta_{DD}^m$ (10.53), and a low $-\log_{10}(\text{pval})$ for the tests related to the flint SNP effect in a flint genetic background. Like for $Vgt2$, a high $-\log_{10}(\text{pval})$ was also detected for tests related to a general SNP effect: $\Delta_{D+F}^m$ (7.47), $\Delta_{DD+DA}^m$ (6.01), $\Delta_{DD+FF}^m$ (7.86) and $\Delta_{DD+DA+FA+FF}^m$ (6.59), but a high $-\log_{10}(\text{pval})$ was detected for the test related to a divergent SNP effect between the dent and the flint genetic backgrounds: $\Delta_{DD-FF}^m$ (3.86). There was also a high $-\log_{10}(\text{pval})$ for a divergent dent SNP effect between different genetic backgrounds: $\Delta_{DD-DA}^m$ (3.03). All these results support the existence of a QTL effect that tends to be higher when the dent genome proportion increases within individuals. It suggests that $Vgt3$ interacts with the genetic background for MF.

The SNP matching a region further referred to as *QTL4.1* on chromosome 4 was detected as associated with MF (20% FDR) using $\Delta_{DD-FF}^m$ (6.93) in $\mathbf{M}_3$. This QTL showed a contrasted effect between alleles of different ancestries with an apparent inversion of effects (Fig. 5-c). This observation was supported by a high $-\log_{10}(\text{pval})$ for the tests related to a divergent SNP effect between ancestries: $\Delta_{D-F}^m$ (5.76), $\Delta_{DD-FF}^m$ (3.48).
<figure><img src="image_6.png" /><figcaption>**Fig 5.** Boxplots of phenotypes adjusted for polygenic background variation using relatedness (MF K corrected) for the different alleles of the six highlighted QTLs: (a) *Vgt2*, (b) *Vgt3*, (c) *QTL4.1*, (d) *QTL2.1*, (e) *QTL7.2* and (f) *Vgt1* using *M₃*. The denomination of the allelic states on the x-axis includes the SNP allele (0/1), its ancestry (D/F) and the genetic background in which it was observed (D/A/F), as presented in Table 2</figcaption></figure>

(6.93) and $\Delta_{(DD+DA)-(FA+FF)}^m$ (5.39). Conversely a low -log₁₀(pval) was detected for tests $\Delta_{DD-DA}^m$ and $\Delta_{FA-FF}^m$, which would have otherwise suggested an interaction with the genetic background. These results support the existence of a local genomic difference at *QTL4.1* between the dent and the flint genetic groups for MF, but no interaction with the genetic background.

The SNP matching a region further referred to as *QTL2.1* on chromosome 2 was detected as associated with MF (5% FDR) using $\Delta_{FA}^m$ (8.24) in *M₃*. This QTL showed a flint effect in the admixed genetic background (Fig. 5-d), which was supported by a high -log₁₀(pval) for the test $\Delta_{FA}^m$ (8.24). Although there was a high -log₁₀(pval) for a general flint SNP effect across genetic backgrounds: $\Delta_{FA+FF}^m$ (5.91), a high -log₁₀(pval) was observed for a divergent SNP effect between those same alleles: $\Delta_{FA-FF}^m$ (3.70). A high -log₁₀(pval) was also observed for a divergent SNP effect between different ancestries in the admixed genetic background: $\Delta_{DA-FA}^m$ (4.84). All these results support the existence of a QTL effect existing only for alleles of flint ancestry in the

352

353

354

355

356

357

358

359

360

361

362

363

364

365
<table><caption>**Table 6.** Information regarding the six highlighted QTLs</caption><thead><tr><th></th><th>Vgt2</th><th>Vgt3</th><th>QTL4.1</th><th>QTL2.1</th><th>QTL7.2</th><th>Vgt1</th></tr></thead><tbody><tr><td>Trait</td><td>MF</td><td>MF</td><td>MF</td><td>MF</td><td>MF</td><td>MF</td></tr><tr><td>SNP</td><td>AX-91100620</td><td>AX-91583310</td><td>AX-91218190</td><td>AX-90601996</td><td>AX-91744673</td><td>MITE</td></tr><tr><td>Chromosome</td><td>8</td><td>3</td><td>4</td><td>2</td><td>7</td><td>8</td></tr><tr><td>Position (Mbp)</td><td>123.50</td><td>158.97</td><td>31.10</td><td>7.04</td><td>173.73</td><td>131.99</td></tr><tr><td>Allele frequency</td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>0DD</td><td>230</td><td>97</td><td>115</td><td>75</td><td>243</td><td>151</td></tr><tr><td>1DD</td><td>70</td><td>203</td><td>185</td><td>225</td><td>57</td><td>149</td></tr><tr><td>0DA</td><td>119</td><td>48</td><td>53</td><td>50</td><td>161</td><td>70</td></tr><tr><td>1DA</td><td>58</td><td>141</td><td>127</td><td>134</td><td>30</td><td>108</td></tr><tr><td>0FA</td><td>81</td><td>92</td><td>107</td><td>74</td><td>113</td><td>17</td></tr><tr><td>1FA</td><td>108</td><td>85</td><td>79</td><td>108</td><td>62</td><td>171</td></tr><tr><td>0FF</td><td>162</td><td>158</td><td>161</td><td>102</td><td>210</td><td>49</td></tr><tr><td>1FF</td><td>142</td><td>146</td><td>143</td><td>202</td><td>94</td><td>255</td></tr><tr><td>-log<sub>10</sub>(pval)</td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>M<sub>1</sub></td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Δ<sup>m</sup> (Dent)</td><td>4.26 *</td><td>10.99 ***</td><td>4.96 *</td><td>0.05</td><td>1.00</td><td>3.34 *</td></tr><tr><td>Δ<sup>m</sup> (Flint)</td><td>2.74 .</td><td>0.88</td><td>0.31</td><td>1.24</td><td>1.20</td><td>0.86</td></tr><tr><td>M<sub>2</sub></td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Δ<sup>m</sup><sub>D</sub></td><td>4.16 *</td><td>9.65 ***</td><td>4.01 *</td><td>0.03</td><td>0.96</td><td>3.37 *</td></tr><tr><td>Δ<sup>m</sup><sub>F</sub></td><td>1.96</td><td>1.10</td><td>2.16 .</td><td>1.29</td><td>0.77</td><td>0.80</td></tr><tr><td>Δ<sup>m</sup><sub>D+F</sub></td><td>5.25 **</td><td>7.47 ***</td><td>0.56</td><td>0.69</td><td>0.10</td><td>0.42</td></tr><tr><td>Δ<sup>m</sup><sub>D-F</sub></td><td>0.57</td><td>3.17 *</td><td>5.76 **</td><td>0.82</td><td>1.47</td><td>3.04 *</td></tr><tr><td>M<sub>3</sub></td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><td>Δ<sup>m</sup><sub>DD</sub></td><td>5.21 **</td><td>10.53 ***</td><td>4.42 *</td><td>0.00</td><td>1.79</td><td>4.62 *</td></tr><tr><td>Δ<sup>m</sup><sub>DA</sub></td><td>2.95 .</td><td>1.38</td><td>1.47</td><td>0.31</td><td>2.64 .</td><td>2.96 .</td></tr><tr><td>Δ<sup>m</sup><sub>FA</sub></td><td>1.09</td><td>2.12 .</td><td>0.97</td><td>8.24 ***</td><td>0.15</td><td>0.15</td></tr><tr><td>Δ<sup>m</sup><sub>FF</sub></td><td>2.85 .</td><td>0.92</td><td>2.34 .</td><td>1.23</td><td>1.51</td><td>0.41</td></tr><tr><td>Δ<sup>m</sup><sub>FA+FF</sub></td><td>2.38 .</td><td>2.00 .</td><td>2.00 .</td><td>5.91 **</td><td>0.49</td><td>0.33</td></tr><tr><td>Δ<sup>m</sup><sub>DD+DA</sub></td><td>5.35 **</td><td>6.01 **</td><td>3.44 *</td><td>0.19</td><td>0.32</td><td>4.96 *</td></tr><tr><td>Δ<sup>m</sup><sub>DA+FA</sub></td><td>3.20 *</td><td>2.93 .</td><td>0.23</td><td>3.09 *</td><td>1.47</td><td>0.85</td></tr><tr><td>Δ<sup>m</sup><sub>DD+FF</sub></td><td>7.15 ***</td><td>7.86 ***</td><td>0.42</td><td>0.70</td><td>0.14</td><td>1.07</td></tr><tr><td>Δ<sup>m</sup><sub>DD+DA+FA+FF</sub></td><td>6.46 **</td><td>6.59 **</td><td>0.39</td><td>2.45 .</td><td>0.63</td><td>1.25</td></tr><tr><td>Δ<sup>m</sup><sub>DA-FA</sub></td><td>0.69</td><td>0.11</td><td>2.11 .</td><td>4.84 *</td><td>2.00 .</td><td>1.52</td></tr><tr><td>Δ<sup>m</sup><sub>DD-DA<sup>b</sup></sub></td><td>0.35</td><td>3.03 *</td><td>0.59</td><td>0.29</td><td>5.58 **</td><td>0.07</td></tr><tr><td>Δ<sup>m</sup><sub>FA-FF</sub></td><td>0.60</td><td>0.69</td><td>0.48</td><td>3.70 *</td><td>1.51</td><td>0.10</td></tr><tr><td>Δ<sup>m</sup><sub>DD-FF</sub></td><td>0.58</td><td>3.86 *</td><td>6.93 **</td><td>0.73</td><td>2.93 .</td><td>2.93 .</td></tr><tr><td>Δ<sup>m</sup><sub>(DD+DA)-(FA+FF)</sub></td><td>0.82</td><td>1.25</td><td>5.39 **</td><td>3.51 *</td><td>0.04</td><td>2.91 .</td></tr><tr><td>Δ<sup>m</sup><sub>(DD+FF)-(DA+FA)<sup>b</sup></sub></td><td>0.73</td><td>0.94</td><td>0.06</td><td>1.35</td><td>1.66</td><td>0.04</td></tr><tr><td>Δ<sup>m</sup><sub>(DD-DA)-(FF-FA)</sub></td><td>0.08</td><td>2.96 .</td><td>0.86</td><td>2.60 .</td><td>6.21 **</td><td>0.13</td></tr></tbody></table>

admixed genetic background. It suggests that *QTL2.1* is specific of flint ancestry and interacts with the genetic background for MF.

The SNP matching a region further referred to as *QTL7.2* on chromosome 7 was detected as associated with MF (20% FDR) using $\Delta_{(DD-DA)-(FF-FA)}^m$ (6.21) in *M3*. This QTL showed contrasted dent effects between the dent and the admixed genetic backgrounds (Fig. 5-e). This observation was supported by a high -log<sub>10</sub>(pval) for the test related to a divergent dent SNP effect between genetic backgrounds: $\Delta_{DD-DA}^m$ (5.58). A high -log<sub>10</sub>(pval) was also observed for the hypothesis testing the equality between the divergent dent SNP effect and the divergent flint SNP effect:

$\Delta_{(DD-DA)-(FF-FA)}^m$ (6.21). All these results support the existence of a QTL with

<footnote>\*\*\* : -log<sub>10</sub>(pval) > 7 ; \*\* : 7 > -log<sub>10</sub>(pval) > 5 ; \* : 5 > -log<sub>10</sub>(pval) > 3 ; . : 3 > -log<sub>10</sub>(pval) > 2</footnote>

<footnote><sup>b</sup> hypothesis testing an interaction between the QTL and the genetic background</footnote>
opposite effects between the dent and the admixed genetic backgrounds. It suggests that *QTL7.2* interacts with the genetic background for MF.

The MITE known to be associated with *Vgt1* was never detected for MF using a FDR of 5% or 20%. However, it showed a dent effect that was conserved between the dent and the admixed genetic background, and no flint effect (Fig. 5-f). This observation is supported by a high -log₁₀(pval) for tests related to the dent SNP effect: Δ^m(Dent) (3.34), Δ^m_D (3.37), Δ^m_D̄ (4.62) and Δ^m_D̄+D̄ (4.96), and a low -log₁₀(pval) for tests related to flint SNP effects. These results support the existence of a local genomic difference at *Vgt1* between flint and dent genetic groups but no interaction with the genetic background for MF.

## Discussion

The stratification of the population sample into distinct genetic groups is a common feature in GWAS studies. Such structure challenges the methods to detect QTLs because (i) spurious associations may be detected if the genetic structure is not accounted for by the statistical model, (ii) QTLs whose polymorphism is correlated with the genetic structure generally have a low probability of being detected when structure/relatedness is accounted for in the model, and (iii) group differences in LD, group-specific genetic mutations and/or epistatic interactions with the genetic background may prevent the detection of SNPs when testing only their average effect.

## Accounting for genetic groups in GWAS

A simple way to deal with genetic groups is to analyze them separately. In our study, a standard GWAS model $M_1$ was applied separately within the dent and the flint datasets. High heritabilities were estimated for each genetic group in the phenotypic analysis, highlighting the suitability of these datasets to detect QTLs. Among the QTLs detected for MF, only one was detected in both dent and flint datasets, and not at the same SNPs, while none were detected in common for FF. One may question whether observing such differences between datasets indicated group specific allele effects, or simply group differences in terms of statistical power. This question often arises when GWAS is applied separately to genetic groups, as in maize [15, 55] or dairy cattle [56, 57], and is very difficult to answer except for obvious configurations such as associations at SNPs segregating only in one group.

Another way to handle genetic groups is to analyze them jointly. One possibility is to apply model $M_1$ while specifying genetic structure as a global fixed effect, in order to prevent the detection of spurious associations. In dairy cattle, this strategy generally improved the precision concerning QTL locations by taking advantage of the low LD extent observed in multi-group datasets. However, while [33] and [32] observed a gain in statistical power due to a larger population size, [31] detected less QTLs by combining breeds compared to separate analyses. They attributed this finding to the limited amount of QTLs segregating within both Holstein and Jersey breeds, but also reported that QTLs detected in both breeds showed only small to medium correlations between within-breed estimates of SNP effects (e.g. 0.082 for milk yield). Obviously, applying $M_1$ jointly to genetic groups does not address directly the problem of whether QTL effects are conserved or not between genetic groups.

A model specifying group specific allele effects was referred to as $M_2$ in this study. As with $M_1$, the existence of a dent ($\Delta^m_D$) and a flint ($\Delta^m_F$) SNP effects can be tested, but $M_2$ also allows us to test the existence of a general ($\Delta^m_{D+F}$) and a divergent ($\Delta^m_{D-F}$) SNP effects between flint and dent ancestries. Note that testing $\Delta^m_{D+F}$ is similar, although not strictly equivalent, to testing a SNP effect by applying $M_1$ to a
multi-group dataset. Using the hypotheses specifically tested in $M_2$ ($\Delta_{D+F}^m$ and $\Delta_{D-F}^m$), it was possible to detect new QTLs that were not detected with $M_1$. In particular, QTLs were detected as having a divergent SNP effect between the dent and flint genetic groups, proving the existence of group-specific QTL effects in this dataset. Several QTLs were detected in common with $M_1$ but each strategy allowed the detection of specific QTLs, demonstrating the complementarity between the models. For equivalent tests in $M_1$ and $M_2$ (e.g. $\Delta^m$ (Dent) in $M_1$ and $\Delta_D^m$ in $M_2$), the lower number of associations detected with $M_2$ can mostly be attributed to a more conservative filtering on allele frequencies. In conclusion, $M_2$ was efficient to identify QTLs with either conserved or specific allele effects between ancestries, but observing group-specific allele effects provided little insight regarding the cause of this specificity. Admixed individuals can help to tackle that issue.

## Benefits from admixed individuals

Admixed individuals were generated for this study by mating pure individuals of each group according to a sparse factorial design. Integrating these admixed individuals in GWAS can be done by simply analyzing the joint multi-group dataset using $M_1$ or $M_2$, which may lead to a gain in statistical power, due to an increase in population size. More interestingly, admixed individuals can be used to disentangle the factors causing the heterogeneity of allele effects across groups.

We developed model $M_3$ to distinguish the allele ancestry (dent/flint) and the genetic background (dent/flint/admixed). 41 QTLs were detected for MF (20% FDR). While many of these QTLs were previously detected using $M_1$ and $M_2$, the new hypotheses tested allowed us to discover new interesting regions. These new QTLs resulted from a gain in statistical power by (i) testing an overall SNP effect for SNP with conserved effects across ancestries and/or genetic backgrounds, or by (ii) testing hypotheses for complex configurations between allele effects. The new hypotheses tested with $M_2$ and $M_3$ did not lead to an increase in false positive rate, based on the observation of the QQ-plots of the test p-values (results not shown).

Note that the idea of exploiting admixed individuals has been proposed in the creation of NAM [39] and MAGIC [40] populations. Compared to our approach, such experimental populations include a limited number of founders, generally selected in different genetic groups. This is beneficial to increase power of detection for alleles which were rare in the initial genetic group(s). However these populations cannot address the question of the epistatic interaction with the genetic background of the original groups, which is possible in our study thanks to the use of numerous parents. Both our approach and NAM and MAGIC designs are therefore expected to have complementary properties.

## Heterogeneity of maize flowering QTL allele effects

From a global perspective, a high number of QTLs have been detected in previous maize studies [15, 21, 36, 58, 59]. When evaluating the American and European NAMs, [21] and [60] showed that flowering time is a trait controlled by a large number of QTLs, many of which display variable effects across individual recombinant populations. Our study highlighted consistently a high number of QTLs and confirmed a large variation in effects. It provides further elements on the origin of this variation, by identifying QTLs affected by local genomic differences, epistasis with the genetic background, or both.

When doing GWAS in a multi-group population, geneticists generally assume that QTL effects are conserved between groups. Such QTLs were detected in our study with the example of the SNP associated with MF in the vicinity of *Vgt2* [14] and its candidate gene: the flowering activator *ZCN8* [61–63] on chromosome 8. At this SNP,
all hypotheses that tested a general SNP effect had a high -log₁₀(pval), and conversely for hypotheses testing a divergent SNP effect. When simultaneously interpreting all tests, *Vgt2* appeared to have an effect that is conserved between genetic groups. Such a QTL can easily be detected in a multi-group population sample using a standard GWAS model [1]. However many QTLs show more complex patterns.

When group-specific allele effects are only due to group differences in LD or group-specific mutations at the QTL, the difference in allele effects should be conserved between the pure and the admixed genetic backgrounds. A first QTL matching this situation (*QTL4.1*) was detected by a SNP located on chromosome 4. High -log₁₀(pval) were observed for the test to a divergent SNP effect between ancestries ($\Delta_m^{(DD+DA)-(FA+FF)}$), suggesting a local genomic difference. To validate this hypothesis, one could produce near isogenic lines with the two alleles from both ancestries introgressed in a dent and a flint genetic backgrounds. A phenotypic evaluation of these individuals would give a definitive proof of a local genomic difference. Nevertheless, it remains difficult to disentangle the effect of LD from that of a genetic mutation without complementary analysis. LD was shown to be different between groups, with a higher LD extent in the dent group (S4 Fig), while LD phases appeared well-conserved at short distances (S5 Fig). However, a strong overall conservation of LD phases at short distances does not exclude a specific configuration for a given SNP-QTL pair. The position of *QTL4.1* is close (< 700 Kbp) to *GRMZM2G126253*, a candidate gene for maize flowering time proposed by [59]. This gene codes for a cullin 3B protein involved in ubiquitination that was shown to be essential to plant development in *Arabidopsis* [64].

Another example is the MITE that we selected based on the *a priori* knowledge that it is associated with *Vgt1* [21, 45, 46] and its candidate gene *ZmRap2.7*. A high -log₁₀(pval) was observed for a dent SNP effect ($\Delta_m^{DD}$ and $\Delta_m^{DD+DA}$) but not for a flint SNP effect. Note that another SNP (AX-91103145) was detected close to the MITE (548 Kbp further), based on 20% a FDR, for $\Delta_m^{(DD+DA)-(FA+FF)}$ (see *QTL8.4* in S7 Fig-*a* and S4 Table). This SNP also showed evidences for a contrasted QTL effect between the dent and flint groups due to a local genomic difference. However these two loci were in very low LD with each other (below 0.05). We can reasonably suggest that the MITE and the SNP both capture a partial but different genetic information of the causal genetic variant at *Vgt1*. [46] already showed the existence of other genetic variants being more associated with maize flowering than the MITE in the vicinity of *Vgt1*, such as CGindel587.

Group-specific allele effects may also be due to an interaction with the genetic background. A first QTL matching this profile was detected by a SNP in the vicinity of *Vgt3* on chromosome 3 [53, 54] and its candidate gene *ZmMADS69* [65]. This QTL showed an effect varying according to the genetic background: large in the dent, intermediate in the admixed and small in the flint. A high -log₁₀(pval) was observed for tests that supported this hypothesis: a dent SNP effect in the dent genetic background ($\Delta_m^{DD}$) and a divergent dent SNP effect between genetic backgrounds ($\Delta_m^{DD-DA}$). If this interaction with the background involves numerous loci, introgressing alleles from a dent to a flint genetic background may lead to disappointing results, as the effect would probably vanish with repeated back-cross generations. If interactions mostly involve a single locus, the effect at *Vgt3* effect is conditioned by the allele at the other locus, so that a simultaneous introgression may be necessary to reach the desired effect. Using near isogenic lines that cumulated an early mutation at *Vgt1* [66] and the early allele at *Vgt3*, the effect of *Vgt3* was shown to vanish in presence of the early allele of *Vgt1* (A. Charcosset pers. comm.), which supports the hypothesis of *Vgt3* interacting with the genetic background. Recently, [65] demonstrated the action of *ZmMADS69*, the candidate gene of *Vgt3*, as being an activator of the regulatory module *ZmRap2.7* -
ZCN8, which are the candidate genes of *Vgt1* and *Vgt2*, respectively. The existence of such interactions is consistent with flowering time being controlled by a network of interacting loci, as now well established in model species arabidopis [67].

Other examples of QTLs interacting with the genetic background were identified. Two of them featured a similar profile in the sense that they mainly exhibited a QTL effect in the admixed genetic background. One was located on chromosome 2 (*QTL2.1*) and showed a flint effect in the admixed genetic background, while the other QTL was located on chromosome 7 (*QTL7.2*) and showed an opposite dent effect between the dent and the admixed genetic backgrounds. Such QTLs are interesting as they are mainly revealed when creating admixed genetic material. They also suggest complex epistatic interactions between QTLs for these traits. The position of *QTL2.1* is close (< 1.4 Mbp) to *ereb197*, a candidate gene for maize flowering time proposed by [59]. This gene codes for an AP2-EREBP transcription factor: a family of transcription factors known to play a role in plant development and response to environmental stress [68]. The position of *QTL7.2* is close (< 100 Kbp) to *dof47*, a candidate gene for maize flowering time proposed by [59]. This gene codes for a C2C2-Dof transcription factor: a family of transcription factors known to play major roles in plant growth and development [69].

The existence of epistatic interactions was also evaluated globally by a test that aimed at detecting directional epistasis [49]. This test was specifically developed to benefit from our admixed genetic material and revealed important directional epistasis for both flowering traits with admixed lines flowering closer to the dent than the flint group. Such epistasis may imply that (i) the effects of early alleles from flint origin tend to decrease in presence of alleles that are more frequent in dent than in flint group and/or (ii) the effect of late alleles from dent origin tends to be promoted by alleles that are more frequent in flint than in dent group. Alternatively, this epistasis can be interpreted as late QTL alleles (common in dent lines but rare in flint lines) interacting in a duplicate way [70], i.e. the presence of a late allele at one QTL is sufficient to confer a late phenotype. This hypothesis is equivalent to early QTL alleles (common in flint lines but rare in dent lines) interacting in a complementary way [70], i.e. early alleles are needed at both loci to confer an early phenotype. We also tested global epistasis that is not directional by decomposing the genetic variance into an additive and an epistatic component, as suggested by [71]. This confirmed the existence of epistatic interactions for FF and MF (S5 Table). In conclusion, the assessment of global epistasis supported the possibility of QTLs interacting with the genetic background, resulting from epistatic interactions with loci that have differentiated allele frequencies between groups. It would be interesting to test the existence of epistatic interactions between pairs of loci. However, a filtering on crossed allele frequencies between pairs of loci would lead to discard most SNPs from the analysis. Other possibilities would be to apply GWAS procedures that are based on testing the epistatic variance of each SNP against the polygenic background [72-74].

## Conclusion

In this study, we proposed an innovative multi-group GWAS method which accounts and tests for the heterogeneity of QTL allele effects between groups. The addition of admixed individuals to the dataset was useful to disentangle the factors causing the heterogeneity of allele effects, being either a local genomic differences or epistatic interactions with the genetic background. Only homozygous inbred lines were considered in this study, but the method may be easily generalized to heterozygous individuals. Recently many studies focused on the problem of genomic prediction across genetic groups [41, 75-78]. In such scenarios, the stability of QTL effects across genetic
backgrounds is an important factor impacting the prediction accuracy. It is also an
important factor of the relevancy of any marker based diagnostic in complex/structured
populations. Our approach opens new perspectives to investigate this stability in a wide
range of species.

575
576
577
578

Supporting information

**S1 Fig.** **Imputation diagram of admixed lines.** Diagram illustrating the procedure applied to impute admixed DH lines from 15K to 600K SNPs using the parental origin of alleles.

579
580
581
582

S2 Fig. Histogram of dent genome proportion among admixed lines.

**S3 Fig.** Genome-wide selection biases among admixed lines. Absolute difference between observed allele frequency of the reference allele $f_o$ estimated on the admixed lines and their expected value $f_e$ along each chromosome ($|f_o - f_e|$). The expected allele frequencies were computed as the mean of flint and dent allele frequencies estimated on the parental lines by taking into account the contribution of each parent. A cubic smoothing spline was adjusted using the R function "smooth.spline", and plotted in red.

**S4 Fig.** LD extent. LD extent estimated, with a sliding window of physical distances between two pairs of loci, in dent and flint genetic groups using the average of (a) the standard $r^2$ or (b) the $r_K^2$ accounting for relatedness between individuals. A cubic smooth spline was adjusted for each group, using the R function "smooth.spline".

**S5 Fig.** Conservation of LD phases. Conservation of LD phases estimated, with a sliding window of physical distances between pairs of loci, using the correlation (a) between the r of dent and flint groups (or the r_K accounting for relatedness between individuals), and (b) between the signs of r in the dent and flint groups (or the signs of the r_K). A cubic smooth spline was adjusted for method, using the R function "smooth.spline".

**S6 Fig.** Position of QTLs detected for FF. Position of QTLs detected for FF with a FDR of 20% using (b) M₁, (b) M₂ and (c) M₃. The size of the grey dots is proportional to the -log₁₀(pval) of the test at the most significant SNP of the region.

**S7 Fig.** Boxplots of phenotypes adjusted for polygenic background variation using relatedness (MF K corrected) for the different alleles of the six other highlighted QTLs: (a) QTL8.4, (b) QTL10.1, (c) QTL3.5, (d) QTL6.3, (e) QTL8.6 and (f) QTL2.2 using M₃. The denomination of the allelic states on the x-axis include the SNP allele (0/1), its ancestry (D/F) and the genetic background in which it is observed (D/A/F), as presented in Table 2.

**S1 Table.** Parameters estimated in the phenotypic analysis. The lines "Row-Column" refer to the modeling of row and columns as defined by the experimental design. AR1 refers to the autoregressive model AR1, while IID refers to the modeling of row and column as being independent and identically distributed among rows and among columns for a given trial. For more information, see the ASReml-R reference manual by [48]. The mean of each trial *j* (with *j* ∈ {2015, 2016}) was computed following: $\mu_j = \mu + \beta_j + \sum_{k=1}^{3} \frac{N_k}{N} \alpha_k$ where $N_k$ is the number of individuals (genotypes)
in genetic background $k$ (with $k \in \{D, A, F\}$) and $N$ is the total number of individuals. The mean of each genetic background was computed following:

$\mu_k = \mu + \alpha_k + \frac{1}{2} \sum_{j=1}^{2} \beta_j$. The genetic variance $\sigma_{G_k}^2$ of each genetic background $k$ and the GxE variance of each genetic background $k$ in each trial $j$ were also reported. The heritabilities of each genetic background $k$ were computed as:

$$h_k^2 = \sigma_{G_k}^2 \left( \sigma_{G_k}^2 + \frac{1}{4} \sum_{j=1}^{2} \sigma_{(G \times \beta)_{jk}}^2 + \frac{1}{4} \sum_{j=1}^{2} \frac{1}{\bar{r}_j} \sigma_{E_j}^2 \right)^{-1}$$ where $\bar{r}_j$ is the mean number of genotype replicates in trial $j$

**S2 Table. Information regarding significant SNPs for MF.** Information regarding significant SNPs for MF using all GWAS strategies: the name of the SNP, the chromosome on which it is located, its position in bp along the chromosome, the frequency of the allelic state observed in the dataset in which it was tested, the GWAS model applied, the hypothesis tested, the $-\log_{10}(\text{pval})$ of the test and the FDR for which it was declared significant.

**S3 Table. Information regarding significant SNPs for FF.** Information regarding significant SNPs for FF using all GWAS strategies: the name of the SNP, the chromosome on which it is located, its position in bp along the chromosome, the frequency of the allelic state observed in the dataset in which it was tested, the GWAS model applied, the hypothesis tested, the $-\log_{10}(\text{pval})$ of the test and the FDR for which it was declared significant.

**S4 Table. Information regarding the six other highlighted QTLs:** QTL8.4, QTL10.1, QTL3.5, QTL6.3, QTL8.6 and QTL2.2.

**S5 Table. Additive, epistatic and residual variance components for each trait with the p-value (pval) of the epistatic component using a likelihood-ratio LR test.** The existence of epistasis can be investigated using a test based on variance components. The epistatic variance component between pairs of loci was estimated on the joint dent, flint and admixed dataset using the model:

$$Y_l = \mu + G_l + (G \times G)_l + E_l \text{ where } (G \times G)_l \text{ is the global epistatic deviation of line } l, \text{ all others terms being identical to those described in } \mathbf{M}_1 \text{ (Eq. 1). Noting}$$

$$g_e^T = ((G \times G)_1, \dots, (G \times G)_N), \text{ one assumes } g_e \sim N(0, K \circ K \sigma_{(G \times G)}^2) \text{ where } K \circ K \text{ is the Hadamard product of the kinship matrix (Eq. 2) with itself and } \sigma_{(G \times G)}^2 \text{ is the epistatic genetic variance between pairs of loci. This model can be seen as a simplified version of the one proposed by [71], as purely homozygous lines were used. The epistatic variance component was tested using a LR test between this model and the same model without the term } (G \times G)_l.$$

**S1 Appendix. Effect of directional epistasis on the mean of an admixed progeny.**

**S2 Appendix. Interpretation of the test $\Delta_{DD-DA}^m$.** This appendix shows that $\Delta_{DD-DA}^m$ tests for an epistatic interaction between the SNP and the genetic background

## Acknowledgments

We thank Stéphane Nicolas (GQE - Le Moulon) for his contribution to genotypic data assembly. We thank Bernard Lagardère, Jean-René Loustalot (INRA Saint-Martin de Hinx) for their contribution to the panel assembly and the coordination of seed
production, all the breeding companies partners of the Amaizing project for the production of admixed lines and the company Limagrain for the genotyping of admixed lines.

659

660

661

## Funding

662

This research was supported by the "Investissement d'Avenir" project (Amaizing, ANR-10-BTBR-03). SR is jointly funded by the program AdmixSel of the INRA metaprogram SelGen and by the partners of the Amaizing project: Arvalis, Caussade-Semences, Euralis, KWS, Limagrain, Maisadour, RAGT and Syngenta. The last seven partners produced admixed lines and Limagrain genotyped these. Funders approved the publication. Otherwise, the funders had no role in study design, data collection and analysis, or preparation of the manuscript.

663

664

665

666

667

668

669

## Data availability statement

670

Genotypes, Allele ancestries, phenotypes, the R script to run analyses, and summary GWAS statistics will be available from the DRYAD (https://datadryad.org/) database (accession number, doi:). All other relevant data are within the paper and its Supporting Information files.

671

672

673

674

## Competing interests

675

The authors have declared that no competing interests exist.

676

## Author Contributions

677

**Conceptualization:** Rio Simon, Mary-Huard Tristan, Moreau Laurence and Charcosset Alain

678

679

**Data Curation:** Rio Simon, Mary-Huard Tristan, Moreau Laurence, Bauland Cyril, Palaffre Carinne, Madur Delphine, Combes Valérie, Charcosset Alain

680

681

**Formal Analysis:** Rio Simon, Mary-Huard Tristan, Moreau Laurence, Charcosset Alain

682

683

**Funding Acquisition:** Charcosset Alain

684

**Investigation:** Rio Simon, Mary-Huard Tristan, Moreau Laurence, Bauland Cyril, Palaffre Carinne, Madur Delphine, Combes Valérie, Charcosset Alain

685

686

**Methodology:** Rio Simon, Mary-Huard Tristan, Moreau Laurence, Charcosset Alain

687

**Project Administration:** Charcosset Alain

688

**Resources:** Rio Simon, Mary-Huard Tristan, Moreau Laurence, Bauland Cyril, Palaffre Carinne, Madur Delphine, Combes Valérie, Charcosset Alain

689

690

**Software:** Rio Simon, Mary-Huard Tristan, Moreau Laurence, Charcosset Alain

691

**Supervision:** Charcosset Alain

692

**Validation:** Rio Simon, Mary-Huard Tristan, Moreau Laurence, Charcosset Alain

693

**Visualization:** Rio Simon, Mary-Huard Tristan, Moreau Laurence, Charcosset Alain

694

**Writing – Original Draft Preparation:** Rio Simon, Mary-Huard Tristan, Moreau

695
## References

1. Yu J, Pressoir G, Briggs WH, Vroh Bi I, Yamasaki M, Doebley JF, et al. A unified mixed-model method for association mapping that accounts for multiple levels of relatedness. Nature Genetics. 2006;38:203-208. doi:10.1038/ng1702.

2. Price AL, Patterson NJ, Plenge RM, Weinblatt ME, Shadick NA, Reich D. Principal components analysis corrects for stratification in genome-wide association studies. Nature Genetics. 2006;38:904-909. doi:10.1038/ng1847.

3. Rincent R, Moreau L, Monod H, Kuhn E, Melchinger AE, Malvar RA, et al. Recovering Power in Association Mapping Panels with Variable Levels of Linkage Disequilibrium. Genetics. 2014;197(1):375-387. doi:10.1534/genetics.113.159731.

4. Pritchard JK, Przeworski M. Linkage Disequilibrium in Humans: Models and Data. The American Journal of Human Genetics. 2001;69(1):1-14. doi:10.1086/321275.

5. Rogers AR. How Population Growth Affects Linkage Disequilibrium. Genetics. 2014;197(4):1329-1341. doi:10.1534/genetics.114.166454.

6. Sawyer SL, Mukherjee N, Pakstis AJ, Feuk L, Kidd JR, Brookes AJ, et al. Linkage disequilibrium patterns vary substantially among populations. European Journal Of Human Genetics. 2005;13:677-686. doi:10.1038/sj.ejhg.5201368.

7. Evans DM, Cardon LR. A comparison of linkage disequilibrium patterns and estimated population recombination rates across multiple populations. The American Journal of Human Genetics. 2005;76:681-687. doi:10.1086/429274.

8. de Roos APWM, Hayes BJ, Spelman RJ, Goddard ME. Linkage disequilibrium and persistence of phase in Holstein-Friesian, Jersey and Angus cattle. Genetics. 2008;179:1503-1512. doi:10.1534/genetics.107.084301.

9. Porto-Neto LR, Kijas JW, Reverter A. The extent of linkage disequilibrium in beef cattle breeds using high-density SNP genotypes. Genetics Selection Evolution. 2014;46(1):22. doi:10.1186/1297-9686-46-22.

10. Badke YM, Bates RO, Ernst CW, Schwab C, Steibel JP. Estimation of linkage disequilibrium in four US pig breeds. BMC Genomics. 2012;13(1):24. doi:10.1186/1471-2164-13-24.

11. Hao C, Wang L, Ge H, Dong Y, Zhang X. Genetic Diversity and Linkage Disequilibrium in Chinese Bread Wheat (Triticum aestivum L.) Revealed by SSR Markers. PLOS ONE. 2011;6(2):1-13. doi:10.1371/journal.pone.0017279.

12. Van Inghelandt D, Reif JC, Dhillon BS, Flament P, Melchinger AE. Extent and genome-wide distribution of linkage disequilibrium in commercial maize germplasm. Theoretical and Applied Genetics. 2011;123(1):11-20. doi:10.1007/s00122-011-1562-3.

13. Technow F, Riedelsheimer C, Schrag TA, Melchinger AE. Genomic prediction of hybrid performance in maize with models incorporating dominance and population specific marker effects. Theoretical and Applied Genetics. 2012;125(6):1181-1194. doi:10.1007/s00122-012-1905-8.
14. Bouchet S, Servin B, Bertin P, Madur D, Combes V, Dumas F, et al. Adaptation of Maize to Temperate Climates: Mid-Density Genome-Wide Association Genetics and Diversity Patterns Reveal Key Genomic Regions, with a Major Contribution of the Vgt2 (ZCN8) Locus. PLOS ONE. 2013;8(8):1–17. doi:10.1371/journal.pone.0071377.

15. Rincent R, Nicolas S, Bouchet S, Altmann T, Brunel D, Revilla P, et al. Dent and Flint maize diversity panels reveal important genetic potential for increasing biomass production. Theoretical and Applied Genetics. 2014;127(11):2313–2331. doi:10.1007/s00122-014-2379-7.

16. Stryjecki C, Alyass A, Meyre D. Ethnic and population differences in the genetic predisposition to human obesity. Obesity Reviews. 2018;19(1):62–80. doi:10.1111/obr.12604.

17. Tang H. Confronting ethnicity-specific disease risk. Nature Genetics. 2006;38(1):12–15. doi:10.1038/ng0106-13.

18. Helgadottir A, Manolescu A, Helgason A, Thorleifsson G, Thorsteinsdottir U, Gudbjartsson DF, et al. A variant of the gene encoding leukotriene A4 hydrolase confers ethnicity-specific risk of myocardial infarction. Nature Genetics. 2006;38(1):68–74. doi:10.1038/ng1692.

19. Barroso I, Luan J, Wheeler E, Whittaker P, Wasson J, Zeggini E, et al. Population-Specific Risk of Type 2 Diabetes Conferred by HNF4A P2 Promoter Variants. Diabetes. 2008;57(11):3161–3165. doi:10.2337/db08-0719.

20. Neuman RJ, Wasson J, Atzmon G, Wainstein J, Yerushalmi Y, Cohen J, et al. Gene-Gene Interactions Lead to Higher Risk for Development of Type 2 Diabetes in an Ashkenazi Jewish Population. PLOS ONE. 2010;5(3):1–6. doi:10.1371/journal.pone.0009903.

21. Buckler ES, Holland JB, Bradbury PJ, Acharya CB, Brown PJ, Browne C, et al. The Genetic Architecture of Maize Flowering Time. Science. 2009;325(5941):714–718. doi:10.1126/science.1174276.

22. Durand E, Bouchet S, Bertin P, Ressayre A, Jamin P, Charcosset A, et al. Flowering Time in Maize: Linkage and Epistasis at a Major Effect Locus. Genetics. 2012;190(4):1547–1562. doi:10.1534/genetics.111.136903.

23. Evangelou E, Ioannidis JPA. Meta-analysis methods for genome-wide association studies and beyond. Nature Genetics. 2013;14:379–389. doi:10.1038/nrg3472.

24. Li YR, Keating BJ. Trans-ethnic genome-wide association studies: advantages and challenges of mapping in diverse populations. Genome Medicine. 2014;6(10):91. doi:10.1186/s13073-014-0091-5.

25. Ioannidis JPA, Ntzani EE, Trikalinos TA. 'Racial' differences in genetic effects for complex diseases. Nature Genetics. 2004;36(12):1312–1318. doi:10.1038/ng1474.

26. Marigorta UM, Navarro A. High Trans-ethnic Replicability of GWAS Results Implies Common Causal Variants. PLOS Genetics. 2013;9(6):1–13. doi:10.1371/journal.pgen.1003566.

27. Ntzani EE, Liberopoulos G, Manolio TA, Ioannidis JPA. Consistency of genome-wide associations across major ancestral groups. Human Genetics. 2012;131(7):1057–1071. doi:10.1007/s00439-011-1124-4.
28. Cole JB, VanRaden PM, O'Connell JR, Van Tassell CP, Sonstegard TS, Schnabel RD, et al. Distribution and location of genetic effects for dairy traits. Journal of Dairy Science. 2009;92(6):2931-2946. doi:10.3168/jds.2008-1762.

29. Hayes BJ, Pryce J, Chamberlain AJ, Bowman PJ, Goddard ME. Genetic Architecture of Complex Traits and Accuracy of Genomic Prediction: Coat Colour, Milk-Fat Percentage, and Type in Holstein Cattle as Contrasting Model Traits. PLOS Genetics. 2010;6(9):1-11. doi:10.1371/journal.pgen.1001139.

30. Cole JB, Wiggans GR, Ma L, Sonstegard TS, Lawlor TJ, Crooker BA, et al. Genome-wide association analysis of thirty one production, health, reproduction and body conformation traits in contemporary U.S. Holstein cows. BMC Genomics. 2011;12(1):408. doi:10.1186/1471-2164-12-408.

31. Raven LA, Cocks BG, Hayes BJ. Multibreed genome wide association can improve precision of mapping causative variants underlying milk production in dairy cattle. BMC Genomics. 2014;15(1):62. doi:10.1186/1471-2164-15-62.

32. van den Berg I, Boichard D, Lund MS. Comparing power and precision of within-breed and multibreed genome-wide association studies of production traits using whole-genome sequence data for 5 French and Danish dairy cattle breeds. Journal of Dairy Science. 2016;99(11):8932-8945. doi:10.3168/jds.2016-11073.

33. Sanchez MP, Govignon-Gion A, Croiseau P, Fritz S, Hozé C, Miranda G, et al. Within-breed and multi-breed GWAS on imputed whole-genome sequence variants reveal candidate mutations affecting milk protein composition in dairy cattle. Genetics Selection Evolution. 2017;49(1):68. doi:10.1186/s12711-017-0344-z.

34. Flint-Garcia SA, Thuillet AC, Yu J, Pressoir G, Romero SM, Mitchell SE, et al. Maize association population: a high-resolution platform for quantitative trait locus dissection. The Plant Journal. 2005;44(6):1054-1064. doi:10.1111/j.1365-313X.2005.02591.x.

35. Camus-Kulandaivelu L, Veyrieras JB, Madur D, Combes V, Fourmann M, Barraud S, et al. Maize Adaptation to Temperate Climate: Relationship Between Population Structure and Polymorphism in the Dwarf8 Gene. Genetics. 2006;172(4):2449-2463. doi:10.1534/genetics.105.048603.

36. Romay MC, Millard MJ, Glaubitz JC, Peiffer JA, Swarts KL, Casstevens TM, et al. Comprehensive genotyping of the USA national maize inbred seed bank. Genome Biology. 2013;14(6):R55. doi:10.1186/gb-2013-14-6-r55.

37. Rius M, Darling JA. How important is intraspecific genetic admixture to the success of colonising populations? Trends in Ecology & Evolution. 2014;29(4):233-242. doi:https://doi.org/10.1016/j.tree.2014.02.003.

38. Brandenburg JT, Mary-Huard T, Rigaill G, Hearne SJ, Corti H, Joets J, et al. Independent introductions and admixtures have contributed to adaptation of European maize and its American counterparts. PLOS Genetics. 2017;13(3):1-30. doi:10.1371/journal.pgen.1006666.

39. McMullen MD, Kresovich S, Villeda HS, Bradbury P, Li H, Sun Q, et al. Genetic Properties of the Maize Nested Association Mapping Population. Science. 2009;325(5941):737-740. doi:10.1126/science.1174320.
40. Cavanagh C, Morell M, Mackay I, Powell W. From mutations to MAGIC: resources for gene discovery, validation and delivery in crop plants. Current Opinion in Plant Biology. 2008;11(2):215-221.
doi:https://doi.org/10.1016/j.pbi.2008.01.002.

41. Rio S, Mary-Huard T, Moreau L, Charcosset A. Genomic selection efficiency and a priori estimation of accuracy in a structured dent maize panel. Theoretical and Applied Genetics. 2019;132(1):81-96. doi:10.1007/s00122-018-3196-1.

42. Bordes J, Dumas de Vaulx R, Lapierre A, Pollacsek M. Haplodiploidization of maize (Zea mays L) through induced gynogenesis assisted by glossy markers and its use in breeding. Agronomie. 1997;17:291-297. doi:10.1051/agro:19970504.

43. Unterseer S, Bauer E, Haberer G, Seidel M, Knaak C, Ouzunova M, et al. A powerful tool for genome analysis in maize: development and evaluation of the high density 600 k SNP genotyping array. BMC Genomics. 2014;15(1):823.
doi:10.1186/1471-2164-15-823.

44. Ganal MW, Durstewitz G, Polley A, BÃ©rard A, Buckler ES, Charcosset A, et al. A Large Maize (Zea mays L.) SNP Genotyping Array: Development and Germplasm Genotyping, and Genetic Mapping to Compare with the B73 Reference Genome. PLOS ONE. 2011;6(12):e1-15.
doi:10.1371/journal.pone.0028334.

45. Salvi S, Sponza G, Morgante M, Tomes D, Niu X, Fengler KA, et al. Conserved noncoding genomic sequences associated with a flowering-time quantitative trait locus in maize. Proceedings of the National Academy of Sciences of the United States of America. 2007;104(27):11376-11381. doi:10.1073/pnas.0704145104.

46. Ducrocq S, Madur D, Veyrieras JB, Camus-Kulandaivelu L, Kloiber-Maitz M, Presterl T, et al. Key Impact of Vgt1 on Flowering Time Adaptation in Maize: Evidence From Association Mapping and Ecogeographical Information. Genetics. 2008;178(4):2433-2437. doi:10.1534/genetics.107.084830.

47. Mangin B, Siberchicot A, Nicolas S, Doligez A, This P, Cierco-Ayrolles C. Novel measures of linkage disequilibrium that correct the bias due to population structure and relatedness. Heredity. 2011;108:285-291.

48. Butler DG, Cullis BR, Gilmour AR, Gogel BJ, Thompson R. ASReml-R Reference Manual Version 4; 2009. VSN International Ltd, Hemel Hempstead, HP1 1ES, UK.

49. Le Rouzic A. Estimating directional epistasis. Frontiers in Genetics. 2014;5:198.
doi:10.3389/fgene.2014.00198.

50. VanRaden PM. Efficient Methods to Compute Genomic Predictions. Journal of Dairy Science. 2008;91(11):4414-4423. doi:10.3168/jds.2010-3719.

51. Laporte F, Charcosset A, Mary-Huard T. Efficient ReML inference in Variance Component Mixed Models using Min-Max algorithms. 2019. Forthcoming.

52. Benjamini Y, Hochberg Y. Controlling the False Discovery Rate: A Practical and Powerful Approach to Multiple Testing. Journal of the Royal Statistical Society. 1995;57(1):289-300. doi:10.2307/2346101.

53. Salvi S, Corneti S, Bellotti M, Carraro N, Sanguineti MC, Castelletti S, et al. Genetic dissection of maize phenology using an intraspecific introgression library. BMC plant biology. 2011;11:4. doi:10.1186/1471-2229-11-4.
54. Salvi S, Emanuelli F, Soriano JM, Zamariola L, Giuliani S, Bovina R, et al. Cloning of Vgt3, a major QTL for flowering time in maize. In: 59th Annual Maize Genetics Conference; 2017.

55. Revilla P, Rodríguez VM, Ordás A, Rincent R, Charcosset A, Giauffret C, et al. Association mapping for cold tolerance in two large maize inbred panels. BMC Plant Biology. 2016;16(1):127. doi:10.1186/s12870-016-0816-2.

56. Buitenhuis B, Janss LL, Poulsen NA, Larsen LB, Larsen MK, Sørensen P. Genome-wide association and biological pathway analysis for milk-fat composition in Danish Holstein and Danish Jersey cattle. BMC Genomics. 2014;15(1):1112. doi:10.1186/1471-2164-15-1112.

57. Buitenhuis B, Poulsen NA, Larsen LB, Sehested J. Estimation of genetic parameters and detection of quantitative trait loci for minerals in Danish Holstein and Danish Jersey milk. BMC Genetics. 2015;16(1):52. doi:10.1186/s12863-015-0209-9.

58. Chardon F, Virlon B, Moreau L, Falque M, Joets J, Decousset L, et al. Genetic architecture of flowering time in maize as inferred from quantitative trait loci meta-analysis and synteny conservation with the rice genome. Genetics. 2004;168(4):2169-2185. doi:10.1534/genetics.104.032375.

59. Li Yx, Li C, Bradbury PJ, Liu X, Lu F, Romay CM, et al. Identification of genetic variants associated with maize flowering time using an extremely large multi-genetic background population. The Plant Journal. 2016;86(5):391-402. doi:10.1111/tpj.13174.

60. Giraud H, Lehermeier C, Bauer E, Falque M, Segura V, Baulaud C, et al. Linkage Disequilibrium with Linkage Analysis of Multiline Crosses Reveals Different Multiallelic QTL for Hybrid Performance in the Flint and Dent Heterotic Groups of Maize. Genetics. 2014;198(4):1717-1734. doi:10.1534/genetics.114.169367.

61. Meng X, Muszynski MG, Danilevskaya ON. The FT-Like ZCN8 Gene Functions as a Floral Activator and Is Involved in Photoperiod Sensitivity in Maize. The Plant Cell. 2011;23(3):942-960. doi:10.1105/tpc.110.081406.

62. Lazakis CM, Coneva V, Colasanti J. ZCN8 encodes a potential orthologue of Arabidopsis FT florigen that integrates both endogenous and photoperiod flowering signals in maize Journal of Experimental Botany. 2011;62(14):4833-4842. doi:10.1093/jxb/err129.

63. Guo L, Wang X, Zhao M, Huang C, Li C, Li D, et al. Stepwise cis-Regulatory Changes in ZCN8 Contribute to Maize Flowering-Time Adaptation Current Biology. 2018;28(18):3005-3015. doi:10.1016/j.cub.2018.07.029.

64. Gingerich DJ, Gagne JM, Salter DW, Hellmann H, Estelle M, Ma L, Vierstra R. Cullins 3a and 3b Assemble with Members of the Broad Complex/Tramtrack/Bric-a-Brac (BTB) Protein Family to Form Essential Ubiquitin-Protein Ligases (E3s) in Arabidopsis. Journal of Biochemical Chemistry. 2005;280(19):18810-18821. doi: 10.1074/jbc.M413247200.

65. Liang Y, Liu Q, Wang X, Huang C, Xu G, Hey S, et al. ZmMADS69 functions as a flowering activator through the ZmRap2.7-ZCN8 regulatory module and contributes to maize flowering time adaptation. New Phytologist. 2019;221(4):2335-2347. doi:10.1111/nph.15512.
66. Chardon F, Hourcade D, Combes V, Charcosset A. Mapping of a spontaneous mutation for early flowering time in maize highlights contrasting allelic series at two-linked QTL on chromosome 8. Theoretical and Applied Genetics. 2005;112(1):1-11. doi:10.1007/s00122-005-0050-z.

67. Bouché F, Lobet G, Tocquin P, Périlleux C. FLOR-ID: an interactive database of flowering-time gene networks in Arabidopsis thaliana. Nucleic Acids Research. 2015;44(D1):D1167-D1171. doi:10.1093/nar/gkv1054.

68. Riechmann J, Meyerowitz E. The AP2/EREBP family of plant transcription factors. Biological chemistry. 1998;379(7):633-646. doi:10.1515/bchm.1998.379.6.633.

69. Yanagisawa S. Dof Domain Proteins: Plant-Specific Transcription Factors Associated with Diverse Phenomena Unique to Plants. Plant and Cell Physiology. 2004;45(4):386-391. doi:10.1093/pcp/pch055.

70. Mather K. Complementary and duplicate gene interactions in biometrical genetics. Heredity. 1967;22:97 EP -. doi:10.1038/hdy.1967.8.

71. Vitezica ZG, Legarra A, Toro MA, Varona L. Orthogonal Estimates of Variances for Additive, Dominance, and Epistatic Effects in Populations. Genetics. 2017;206(3):1297-1307. doi:10.1534/genetics.116.199406.

72. Jannink JL. Identifying Quantitative Trait Locus by Genetic Background Interactions in Association Studies. Genetics. 2007;176(1):553-561. doi:10.1534/genetics.106.062992.

73. Crawford L, Zeng P, Mukherjee S, Zhou X. Detecting epistasis with the marginal epistasis test in genetic mapping studies of quantitative traits. PLOS Genetics. 2017;13(7):1-37. doi:10.1371/journal.pgen.1006869.

74. Legarra A, Vitezica ZG, Naval-Sánchez M, Henshall J, Raidan F, Li Y, et al. Association analysis of loci implied in “buffering” epistasis. bioRxiv. 2019;637579. doi:10.1101/637579.

75. de Roos APW, Hayes BJ, Goddard ME. Reliability of Genomic Predictions Across Multiple Populations. Genetics. 2009;183(4):1545-1553. doi:10.1534/genetics.109.104935.

76. Chen L, Schenkel F, Vinsky M, Crews DH, Li C. Accuracy of predicting genomic breeding values for residual feed intake in Angus and Charolais beef cattle. Journal of Animal Science. 2013;91(10):4669-4678. doi:10.2527/jas.2013-5715

77. Guo Z, Tucker DM, Basten CJ, Gandhi H, Ersoz E, Guo B, et al. The impact of population structure on genomic prediction in stratified populations. Theoretical and Applied Genetics. 2014;127(3):749-762. doi:10.1007/s00122-013-2255-x

78. Duhnen A, Gras A, Teyssèdre S, Romestant M, Claustres B, DAYDE J, et al. Genomic Selection for Yield and Seed Protein Content in Soybean: A Study of Breeding Program Data and Assessment of Prediction Accuracy. Crop Science. 2017;57(3):1325-1337. doi:10.2135/cropsci2016.06.0496.