Figure 1. Populus tremula leaf physiology and
morphology. A) A representative mature terminal leaf.B) A representative flushed short shoot bud with mature
pre-formed leaves. C) Cross-section (15 µm thick) of a
representative mature terminal leaf. D) Cross-section (10 µm
thick) of a representative mature pre-formed leaf. Cross sections were
embedded in 1,5 % plant agarose, sectioned using a vibratome and
stained for one minute with Toluidine Blue; Magnification 200X.
Figure 2 . An overview of leaf circularity, indent depth and
leaf area in the Swedish Aspen (SwAsp) collection. A) Example
leaves from the SwAsp collection. Each leaf is from one genotype and the
set of leaves was selected to reflect the extent of variation
represented among all genotypes. B) Genetic correlations among
measures of leaf circularity, median indent depth and leaf area for data
recorded in two common gardens (Ekebo, southern Sweden and Sävar,
northern Sweden) and two years (2008 and 2011). C) Heatmap
representations of leaf circularity (left), median indent depth (centre)
and leaf area (right) in the two common gardens and years. Each
represented value is a median calculated from 10 leaves per clonal
replicate (median n=3) of each genotype. In each heat map, genotypes are
clustered by Euclidean distance with the population from which the
genotype originates indicated by coloured bars to the left of the heat
map. The values represented were used to perform the genome-wide
association mapping results depicted in D. D) Intersection of
the top 1,000 ranked SNPs from genome wide association mapping for the
three traits measured in the two years and gardens.
Figure 3. The genomic context distribution of SNPs within the
top ranked 1,000 genes identified using association mapping for leaf
circularity and indent depth.
Figure 4. Selection of a set of high and low phenotype value
genotypes. Genotypes with high and low phenotypic values of circularity
(A ) and indent depth (B ) were selected on the basis of
the top or bottom quartiles of the Best Linear Unbiased Prediction
(BLUP) values for each trait. The colour intensity reflects the four
quartiles of trait values with darker colours meaning higher values. The
left density plots represent BLUPs for the two traits. The distribution
plots in the centre depict the original trait measurements in the two
gardens and the two years of sampling prior to the BLUP calculation.
Gene expression assayed in young leaf buds (Mähler et al. 2017)
was used to test for differential expression between genotypes within
the low and high sets for each phenotypic trait. The Venn diagram shows
the intersection between the significantly differentially expressed
genes and the top 1,000 genes from the genome wide association study
(GWAS) for the two traits.
Figure 5. Overview of gene expression in developmental series
of terminal leaves. A) Principal component analysis of the
1,000 most variable genes. Data points are shaded by tissue type (Pf =
Pre-formed; T = terminal) and leaf number in addition to different point
styles being used. Ellipses around each group indicate a 95% confidence
area for a particular sample group. B) Differential expression
matrix of all sample points showing all significantly up-regulated genes
(P adj < 0.05) in the lower triangular
matrix and all significantly down-regulated genes
(P adj < 0.05) in the upper triangular
matrix. In the lower triangular, columns constitute the numerator of the
differential test, while rows constitute the numerator in the upper
triangular half. T = terminal.
Figure 6. Visualisation of differential expression within the
terminal developmental series gene co-expression network. A)Transition between stage -4 and -3. B) Transition between stage
-3 and -2. C) Transition between stage -2 and -1. D)Transition between stage -1 and 0. E) Transition between stage
0 and 1. F) Transition between stage 1 and 2. Nodes are sized
inversely by adjusted P -value of a differential expression
comparison at the relevant transition and shaded by the fold change
(higher in red, lower in blue) with the earlier leaf number as the
numerator. P- values ranged from 1.519197e-104 to 9.999962e-01
with the highest P -value nodes being 100x smaller than the lowest
ones.
Figure 7. Developmental expression of genes differentially
expressed between phenotypic extremes. The heatmap shows the expression
of the differentially expressed genes depicted in Figure 4. Expression
in the terminal leaf development data set is represented for genes
differentially expressed between phenotypically extreme genotype sets
for circularity (A ) and indent depth (B ). Genes that
were found among the top 1,000 GWAS genes are indicated by coloured
lines on the left side of the heatmap. Genes are clustered hierarchical
clustering using average linkage. Samples are arranged by leaf age.
Values depicted are VST (Variance-Stabilising Transformation) normalised
expression. Clustering was performed using complete linkage and
correlation-based distance.
Supplementary Figure S1 . Boxplots showing genotype x
environment interactions for size and shape metrics in pre-formed leaves
of the SwAsp collection in common gardens at Ekebo (56 °N) and Sävar (62
°N). Each box represents a common garden in a given year with 111 to 113
genotypes in clonal replicates. (A) Leaf circularity, (B) indent depth,
and (C) leaf area. GxE analysis of variance F and P -values
are provided in Table 1 Boxes annotated with the same letter are not
significantly different at P < 0.05.
Supplementary Figure S2 . Manhattan and QQ-plots for leaf area,
circularity, and indent depth.
Supplementary Figure S3 . Venn diagram of the top 1,000 genes
from the genome wide association analysis for the three leaf traits.
Supplementary Figure S4 . Correlation between gene expression
among genotypes of the Swedish Aspen collection and leaf shape and size
traits. The left panel shows the distribution for all expressed genes
for each phenotypic trait. The right panel separates genes that were
differentially expressed between sets of genotypes with contrasting
trait values from non-differentially expressed genes.
Supplementary Figure S5. Network clusters of genes within the
community expression network inferred for the terminal leaf
developmental series. Clusters were defined using InfoMap. The
corresponding eigengene expression profile is shown for each partition
as well as representative enriched PFAM and GO categories. * indicates
an adjusted P -value for the terms between 0.05 and 0.01, **
indicates 0.01 - 0.001, *** indicates < 0.001. Only terms in
the same range of significance were included for a cluster (if any). The
node size is scaled according to the node’s PageRank.
Supplementary Figure S6 . Gene set enrichment analysis results
for 1,000 top genes from the genome wide association (GWAS) and
differentially expressed gene sets from the GWAS population (DEGs)
within the terminal developmental gene co-expression network.
Supplementary Figure S7 . Heatmap of the scaled and centred gene
expression values from the terminal leaf development series for all
genes that were differentially expressed between genotype extreme groups
for at least one of the traits, indicated by the green marks on the
left.
Supplementary Figure S8 . Venn diagram of the genes
differentially expressed between phenotypic extremes (high and low BLUP
quartiles) of the three leaf traits.
Supplementary Figure S9. Scaled network centrality for genes
that were differentially expressed between sets of genotypes with
contrasting trait values and non-differentially expressed genes within
the Swedish Aspen, terminal and pre-formed leaf development series
co-expression networks. The Swedish Aspen co-expression network is from
Mähler et al . (2017) and uses WGCNA-calculated network centrality
scores. The leaf development network uses InfoMap scores. Both scores
were scaled to uniform scale to enable comparison. Results were compared
to randomised sets of genes to ensure that the observed signal was not a
general pattern. Genes that were also among the top 1,000 genes
identified by genome wide association mapping are indicated as coloured
data points.
Supplementary Figure S10. Gene set enrichment analysis results
for genome wide association and population differentially expressed gene
sets within the Swedish Aspen gene co-expression network from Mähleret al. (2017).
Supplementary Figure S11. Tajima’s D values for the top 1,000
genes identified by genome wide association mapping and for genes that
were differentially expressed between sets of genotypes with contrasting
trait values. Tajima’s D values were obtained from Mähler et al.(2017). Differentially expressed genes that were also among the top
1,000 genes identified by genome wide association mapping are indicated
as coloured data points.
Supplementary Table S1 . Top 1,000 genome wide associations for
leaf area, circularity, and indent depth. The columns are snp: SNP name;
chr: scaffold where SNP is located; pos: base-pair position of SNP on
scaffold; pseudo_chr: pseudo-chromosome where SNP is located;
pseudo_pos: base-pair position of SNP on pseudo-chromosome;
minor_allele: minor allele of SNP; major_allele: major allele of SNP;
maf: minor allele frequency of SNP; beta: effect size of the
association; se: standard error of the beta; pvalue: likelihood ratio
test P -value of the beta being different from zero; gene: gene
that the SNP is located in/near (within 2 kbp), if any; feature: gene
feature that the SNP is located in, if any; at_gene: best BLAST hit inArabidopsis thaliana ; potri_gene: best BLAST hit inPopulus trichocarpa .
Supplementary Table S2 . Top 1,000 genes from the genome wide
associations for leaf area, circularity, and indent depth. The genes
were selected by rank ordering the GWAS results by associationP -value and then picking genes from the top of the list until
1,000 unique genes had been selected.
Supplementary Table S3 . Differential expression results for the
three leaf traits. In addition to the statistics from DESeq2, it has
been indicated whether the gene is also an eGene in the SwAsp
collection, whether the gene is part of the leading-edge subset in the
gene set enrichment analysis for the trait, gene annotations, reciprocal
best DIAMOND hit in A. thaliana , as well as the best DIAMOND hit
in A. thaliana .
Supplementary Table S4 . Node centrality statistics of all genes
in the terminal leaf network. The centrality statistics computed are
PageRank, closeness, betweenness, strength, eigenvector, Katz,
Kleinberg’s authority and Kleinberg’s hub score.
Supplementary Table S5 . Gene ontology and PFAM (protein family)
enrichments for all clusters with more than five members in the terminal
leaf network.
Supplementary Table S6 . Normalised and raw expression values
and results of all pairwise differential expression comparisons
containing mean expression, log2 fold changes, standard
error of the log2 fold changes, the test statistic, the
test P value and the adjusted P value for all genes in all comparisons.