Seventeen of these modules had significant Gene Ontology enrichment with each module consisting of transcript abundance pat- terns either predominately regulated by cis- or trans- eQTL . To determine which IL’s are important for module regulation, the median transcript abundance value of module genes for each IL was calculated and used to identify IL’s with significantly altered module steady-state transcript level. Three modules were present in all mappings of the BH-SNE deter- mined through iterations of DBscan analysis and GO enrichment and were designated as landmark modules . The largest module had a GO enrichment for photosynthesis and related processes, and significant trans-eQTL scattered widely across the genome with no bin or IL identified as the primary regulating region . The second landmark module was enriched for transcript abundance patterns with roles in defense, metabolism, and signaling with the majority of their trans-eQTL mapped to IL6.2 and 6.2.2 . The third module, containers size for raspberries which is enriched for transcript abundance patterns with Cys-type peptidase activity, was predominately composed of genes regulated by cis-eQTL on IL 4.2, 4.3, and 4.3.2 .

A cluster of genes enriched for “peptidase regulation” also emerged from a transcriptome study of leaf development for three tomato species; this cluster was uniquely associated with S. pennellii orthologs at the P5 stage of leaf development, indicating that this species has a unique pattern of gene expression, whichinvolves peptidase regulation , and may be related to leaf maturation and senescence processes .One of the landmark modules from the clustering analysis was enriched for transcript abundance pat- terns related to plant defense . Therefore, we explored the genetic basis of transcriptional changes associated with plant defense. IL6.2 and IL6.2.2, and associated bins 6B and 6C, in particular, influence of the transcriptional responses of genes associated with plant defense and signaling . The genes showing increased steady-state transcript levels in both IL’s com- pared to cv M82, as well as the genes regulated by the corresponding bins, show enrichment of the GO categories response to stress and stimulus, cell death, defense response, and plant-type hypersensitive response . Promoter enrichment analysis for these genes showed enrichment of a W-box promoter motif that is recognized by WRKY transcription factors and influences plant defense response . Both bins, in particular bin 6C, contain genes involved in pathogen, disease, and defense response, such as NBS-LLR resistance genes, WRKY transcription factors, Multidrug resistance genes, Pentatricopeptide repeat containing genes, Chitinase, and Heat Shock Protein coding genes. This transcriptional response in the IL’s is also reflected in the morphology of IL6.2.2; the plants are necrotic and dwarfed .

Previously, a phenotypic study for the chromosome 6 introgression, specifically a 190-kb region in bin 6C, in a pathogen /control experiment was shown to confer hypersensitive response in IL6.2 and 6.2.2 . Taken together, these findings suggest bins 6B and 6C contain master genetic regulators of plant defense response genes, though identification of the causal gene/s that influence many other genes in trans will need further genetic dissection of these bins.Given the striking differences in leaf features be- tween S. pennellii and cv M82 that are manifested in many IL’s , the IL population provides an excellent system for determining the extent of genetic regulation of genes controlling leaf development. Previous phenotypic and QTL analyses identified many IL’s, such as IL4.3, IL8.1.5, IL8.1.1, and IL8.1, harboring loci regulating leaf and plant developmental traits . IL4.3, which harbors loci with the largest contribution to leaf shape and shows larger epidermal cell sizes , exhibited decreased steady-state transcript levels for many genes associated with cell division, such as Cyclin-dependent protein kinase regulator-like protein , Cyclin A-like protein , and F-box/LRR repeat protein 2 SKP2A . In addition, genes showing differences in transcript levels in IL4.3 were enriched for the promoter motifs MSA and the E2F binding site . Genes with decreased transcript levels in IL’s 8.1.5, 8.1.1, and 8.1, also included genes associated with leaf development and morphology, genes encoding WD-40 repeat family protein LEUNIG, Homeobox- Leu zipper protein PROTODERMAL FACTOR 2, and the transcription factor ULTRAPETALA . We further investigated the transcript expression dynamics of a set of literature-curated genes related to leaf development across the IL’s and bins . A number of canonical leaf developmental genes such as SHOOTMERISTEMLESS , GROWTH-REGULATING FACTOR 1 , ARGONAUTE 10 , BELL , LEUNIG , and SAWTOOTH 1 were differentially expressed at the transcript level in more than five IL’s.

At the level of bins, genes involved in leaf development were regulated by eQTL scattered widely across the genome . eQTL-regulation of leaf developmental genes for some of IL’s, such as IL 2.1, 4.3, 5.4, 8.1/8.1.1/8.1.5, and 9.1.2 showing strong leaf phenotypes, is summarized in Supplemental Table S4. We then examined the location of literature-curated leaf developmental genes within the identified modules in the BH-SNE mapping . The highest number of literature-curated leaf developmental genes was located in the photosynthesis module, whereas 19 of these genes were located in the leaf development module , suggesting a relationship between these two modules. Over one-third of the transcript ex- pression patterns in the leaf development module have significant eQTL that map to bins 4D, 8A, and 8B , suggesting that these bins contain important regulators of leaf development. This enrichment of eQTL for specific bins is also consistent with the strong leaf phenotypes for IL’s 4.3, 8.1, 8.1.1, and 8.1.5. Altogether, DE, eQTL, and BH-SNE results indicate that while there is no obvious master regulatory bin for leaf developmental genes, many are under strong genetic regulation by eQTL distributed throughout the genome . This observation underscores the highly polygenic regulation of leaf development as multiple loci, residing in many different chromosomal locations, regulate the expression of key leaf- developmental genes at the transcriptional level.Since photosynthesis GO terms were enriched for the largest module from the clustering analysis and there was a correlation between photosynthesis and leaf developmental modules , we examined the genetic regulation of photo- synthetic genes by specific IL’s and corresponding bins. Genes related to photosynthesis show increased transcript levels across 21 IL’s distributed on all chromosomes except chromosome 5 , showing multigenic regulation of photosynthetic traits. Many of these IL’s, including 8.1.5, 8.1.1, 8.1, and 4.3, and associated bins showed regulation of genes linked to photosynthesis, chlorophyll biosynthesis, and response to light stimulus . This observation indicates that IL’s may also differ from each other and from the cultivated M82 background in photosynthetic efficiency. However, no studies, so far, have investigated the photosynthetic phenotype of these IL’s. To analyze the relationship between the leaf development and photosynthesis modules, the median transcript abundance value of all genes in each module was compared, resulting in a significant negative correlation . This analysis likely reflects the transition from leaf development to leaf maturation captured in our shoot meristem samples. The genes found in the leaf development module may promote developmental processes such as cell division and maintenance or meristematic potential, whereas the leaf development-related genes found in the photosynthesis module may act to suppress this process to allow for maturation of the leaf. The two modules had their most influential eQTL on bins 4D, 8A, and 8B , big plastic pots suggesting that leaf development and photosynthetic genes not only have transcript levels in opposition but also likely share common regulatory loci. This finding is consistent with the link between leaf development and photosynthesis that we established previously by meta- analysis of developmental and metabolic traits .Since the eQTL study used shoot apices that includes the shoot apical meristem and developing leaves, we resolved the detected eQTL to specific tis- sues and temporally regulated development using previous gene expression data. We analyzed transcript abundance in laser micro-dissected samples representing the SAM + P0 versus the P1 that represents transcript levels in the meristem and the first differentiated leaf . We also analyzed hand dissected samples of the SAM + P0-P4 vs. the P5 collected over time , representing genes regulated by vegetative phase change .

Using a bootstrapping approach, we identified bins statistically enriched for genetically regulating genes with previously identified transcript expression patterns . Except for one instance , bins enriched for transcript expression patterns represented trans regulation, hinting at predominant regulation of gene expression patterns mediated by transcription factors at the level of transcription. Most SAM/P0 versus P1 enriched bins were enriched for P1 transcript expression . We previously showed that genes with high P1 transcript levels are enriched for photosynthetic-related GO terms compared to SAM/P0 genes enriched for transcription, cell division, and epigenetics related GO terms , suggesting a genetic basis at both a functional and tissue specific level for genes related to photosynthesis expressed preferentially in the P1 compared to the SAM/P0. Bins enriched for regulation of genes with temporally dependent steady-state transcript levels were mostly associated with genes with decreasing transcript level over time, for both the SAM + P0- P4 and P5 . Interestingly, 3 bins share enrichment for genes with decreasing transcript levels over time in both the SAM + P0-P4 and P5 , suggesting true temporal trans regulation, regardless of tissue, by these loci. Broadly, genes with increasing transcript levels over time are associated with transcription and small RNA GO terms in both the SAM + P0-P4 and P5, whereas decreasing transcript levels over time are associated with translation associated GO terms in the SAM + P0-P4 and photosynthetic activity in the P5 .To connect detected eQTL with leaf and hypocotyl phenotypes under two different environmental conditions, we correlated transcript abundance with leaf number, leaf complexity , and hypocotyl length phenotypes of the IL’s grown under simulated sun and shade conditions. Significant correlations with transcript abundance pat- terns were identified for all three phenotypes analyzed under both treatments . Focusing on a subset of these transcript expression pat- terns that had associated eQTL enabled us to connect the phenotypes to their regulatory loci . Genes negatively correlated with leaf number showed enrichment of leaf development GO terms, whereas positively correlated genes showed enrichment of photosynthesis-related GO terms . For the leaf complexity trait, correlations were reversed compared to leaf number . The transcript levels of these genes associated with leaf number were pre- dominantly regulated by eQTL on chromosomes 7 and 8 and those of leaf complexity on chromosomes 4, 7, and 8 . These results, in combination with DE, eQTL, and BH-SNE, highlight bins on chromosomes 4 and 8 as important genetic regulators of leaf develop- mental genes.Five genes were positively correlated with hypocotyl length under simulated shade, and one gene was negatively correlated with hypocotyl length under both sun and shade . eQTL for the positively correlated genes are located on chromosomes 3, 7, and 11, whereas the single cis-eQTL for the negatively correlated gene, Solyc10g005120 , was located in bin 10A.1 . The transcript is expressed only in IL 10.1, which has the S. pennellii version of the gene and an attenuated shade avoidance response, but is not expressed in IL 10.1.1, which has the M82 version of the gene and a normal shade avoidance response . This indicates that genes in bin 10A, the nonoverlapping regions of 10.1 and 10.1.1, are responsible for the shade avoidance response. Bin 10A includes Solyc10g005120, the one gene negatively correlated with hypocotyl length under both sun and shade. A set of Backcross Inbred Lines , developed from cv M82 and S. pennellii, provide higher resolution gene mapping with smaller bin sizes . To further explore the role of Solyc10g005120, we used BIL-128, which contains a subregion of bin 10A and has a secondary introgression on chromosome 2 . Influence of the secondary introgression was examined using BIL-033, which shares the introgression on chromosome 2. BIL- 128 has an attenuated shade avoidance response, as does 10.1, whereas BIL- 033 undergoes a shade avoidance response similar to that of cv M82 . These results rule out the influence of chromosome 2 genes on the attenuated hypocotyl phenotype and confirm the influence of the bin 10A subregion, which includes Solyc10g005120, on the attenuated hypocotyl phenotype . Solyc10g005120 is an uncharacterized gene, and our observations highlight it as a new candidate regulating shade avoidance responses.Plant materials, growth conditions, and experimental design were described in , but are outlined here briefly. Seeds of wild tomato IL’s and cultivated tomato were obtained either from Dani Zamir or from the Tomato Genetics Resource Center . Seeds were stratified in 50% bleach for 2 min and grown in darkness for 3 d for uniform germination before moving to a growth chambers for 5 d.