Estimation of maternal expression trait values
Prior to conducting further analyses using the initial 47,651
transcripts obtained after a series of standard filtering steps (Methods
S1), we normalised read count data for varying library sizes using the
trimmed mean of M values (TMM) approach as implemented with thecalcNormFactors function from the edgeR v.3.14.0 package
(Robinson et al ., 2010) in R. We further filtered transcripts
that were not expressed across all three siblings nested within the same
maternal tree in either garden. We implemented a linear mixed effects
model for each of the remaining 27,441 transcripts (Eq. 1):
\(Y_{\text{ijklmn}}=\ \mu+\text{Batch}_{j}+\ \text{SY}_{k}+\text{Garden}_{n}|\left(\frac{\text{Pop}_{l}}{\text{Fam}_{m}}\right)+\ \epsilon_{\text{ijklnm}}\)(1)
Here, \(Y_{\text{ijklmn}}\) represents the normalised expression value
of a transcript for the i th seedling in thej th sequencing batch sowed in thek th year from the m thmaternal tree nested within the l th population
at n th garden. The global mean across all
seedings for a transcript is given by \(\mu\). The terms Batchand SY are fixed effects and represent the date of sequencing and
year of seed sowing, respectively. The term Fam represents the
maternal tree that is nested within Pop , which represents the
source population of the seeds for a maternal tree. Both Fam andPop are treated as random effects that are allowed to vary by the
fixed effect of Garden . This linear mixed effect model was fitted
using the fitVarPartModel function with log-transformed
normalised counts and precision weights, which are ideal for RNAseq
datasets (Law et al ., 2014), obtained from thevoomWithDreamWeights function in the variancePartition v.1.16.1
package (Hoffman & Schadt, 2016) in R. Using the fitted model in
equation (1), we obtained garden-specific genetic values for maternal
trees and populations for each transcript.