Get 20M+ Full-Text Papers For Less Than $1.50/day. Start a 7-Day Trial for You or Your Team.

Learn More →

Systems mapping: how to improve the genetic mapping of complex traits through design principles of biological systems

Systems mapping: how to improve the genetic mapping of complex traits through design principles... Background: Every phenotypic trait can be viewed as a “system” in which a group of interconnected components function synergistically to yield a unified whole. Once a system’s components and their interactions have been delineated according to biological principles, we can manipulate and engineer functionally relevant components to produce a desirable system phenotype. Results: We describe a conceptual framework for mapping quantitative trait loci (QTLs) that control complex traits by treating trait formation as a dynamic system. This framework, called systems mapping, incorporates a system of differential equations that quantifies how alterations of different components lead to the global change of trait development and function through genes, and provides a quantitative and testable platform for assessing the interplay between gene action and development. We applied systems mapping to analyze biomass growth data in a mapping population of soybeans and identified specific loci that are responsible for the dynamics of biomass partitioning to leaves, stem, and roots. Conclusions: We show that systems mapping implemented by design principles of biological systems is quite versatile for deciphering the genetic machineries for size-shape, structural-functional, sink-source and pleiotropic relationships underlying plant physiology and development. Systems mapping should enable geneticists to shed light on the genetic complexity of any biological system in plants and other organisms and predict its physiological and pathological states. Background success of this prediction depends on how well we can Predicting the phenotype from the genotype of complex map the underlying QTLs and characterize complex organisms is one of the most important and challenging interactions of these QTLs with each other and with questions we face in modern biology and medicine [1]. environmental factors. Powerful statistical models have Genetic mapping, dissecting a phenotypic trait to its been developed in the past two decades to detect QTLs underlying quantitative trait loci (QTLs) through the and study their biological function in a diverse array of use of molecular markers, has proven powerful for phenotypic traits [3-9]. Worldwide, a substantial effort establishing genotype-phenotype relationships and pre- has been made resulting in the collection of a large dicting phenotypes of individual organisms based on amount of data aimed at the identification of QTLs [10-15]. Unfortunately, despite hundreds of thousands of their QTL genotypes responsible for the trait [2]. The QTLs detected in a diversity of organisms, only a few of them have been isolated by positional cloning (see * Correspondence: [email protected]; [email protected] [16-18]), leaving it unsolved how to construct a geno- † Contributed equally Center for Computational Biology, National Engineering Laboratory for Tree type-phenotype relationship map using genetic mapping. Breeding, Key Laboratory of Genetics and Breeding in Forest Trees and The most likely reason for this result may arise from a Ornamental Plants, Beijing Forestry University, Beijing 100083, China possibility that the QTLs detected by stringent statistical National Center for Soybean Improvement, National Key Laboratory for Crop Genetics and Germplasm Enhancement, Soybean Research Institute, tests are not biologically relevant. Existing strategies for Nanjing Agricultural University, Nanjing 210095, China QTL mapping were built on testing for a direct Full list of author information is available at the end of the article © 2011 Wu et al; licensee BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited. Wu et al. BMC Systems Biology 2011, 5:84 Page 2 of 11 http://www.biomedcentral.com/1752-0509/5/84 association between genotype and end-point phenotype. production as well as internal efficiency through redu- Although such strategies are simple and have been cing distances and the time to transport water, nutri- widely accepted, they neglect the biological processes ents, and carbon. The integration of this optimization involved in trait development [14]. To attempt to fill theory expressed in terms of allometric scaling with the this gap, a statistical model, called functional mapping coordination theory leads to a tripled group of ordinary [19-21], has been developed to study the interplay differential equations (ODEs) to specify the coordination of the leaf, stem, and root biomass for a plant: between genetics and the developmental process of a phenotypic trait by integrating mathematical models and dM computational algorithms. If a trait is understood as a = α W − γ M L L L dt “system” that is composed of many underlying biological dM components [22-24], we should be in a better position β (1) = α W dt to comprehend the process and behavior of trait forma- dM tion based on interactive relationships among different R = α W − γ M R R R components. Through mapping and using those QTLs dt that govern design principles of a biological system, a new trait that is able to maximize resource-use effi- where M , M ,and M are the biomasses of the leaves L S R ciency can be generated and engineered. (L), the stem (S), and the roots (R), respectively, with As one important strategy for plants to respond to whole-plant biomass W = M + M + M ; a and b are L S R variation in the availability of resources in their environ- the constant and exponent power of an organ biomass ment, biomass allocation has been extensively used to scaling as whole-plant biomass [32,33]; and g is the rate study the relationship between structure and function in of eliminating ageing leaves and roots. The complex modern ecology [25-28]. The concept of biomass alloca- interactions between different parts of a plant that tion has now been increasingly integrated with plant underlie design principles of plant biomass growth can management and breeding, aimed to direct a maximum be modeled and studied by estimating and testing the amount of biomass to the target of harvest (leaves, stem, ODE parameters (a , b , l , a , b , a , b , l ). For L L L S S R R R roots, or fruits) [29-31]. If the whole-plant biomass is example, plants are equipped with a capacity to optimize considered as a target trait, we need to understand how their fitness under low nutrient availability by shifting different organs of a plant coordinate and interact to the partitioning of carbohydrates to processes associated optimize the capture of nutrients, light, water, and car- with nutrient uptake at a cost of carbon acquisition bon dioxide in a manner that maximizes plant growth [29]. These parameters can be used to quantify and pre- rate through a specific developmental program because dict such regulation between different plant parts in plant biomass growth is not simply the addition of bio- response to environmental and developmental changes. mass to various organs. Many theories and models have In this article, we put forward a conceptual framework been proposed to predict the pattern of biomass parti- to incorporate the design principles of trait formation tioning in a response to changing environment. Chen and development into a statistical framework for QTL and Reynolds [27] used coordination theory to model mapping. Complementary to our previous functional the dynamic allocation of carbon to different organs mapping [19-21], we name this new mapping framework during growth in relation to carbon and water/nitrogen “systems mapping” in light of its systems dissection and supply by a group of differential equations. Compared modeling of phenotypic formation. A group of ODEs to the conventional optimization model in the context like (1) or other types of differential equations is used to of maximizing the relative growth rate of a plant, the quantify the phenotypic system. Much work in solving coordination model does not require an unrealistic ODEs has focused on the simulation and analysis of the capacity the plant possesses for knowing beforehand the behavior of state variables for a dynamic system, but the environmental conditions it will experience during the estimation of ODE parameters that define the system growth period. Here, we integrate the coordination and based on the measurement of state variables at multiple optimization model to study the pattern of biomass par- time points is relatively a new area. Yet, in the recent titioning by incorporating the allometric scaling theory years, many statisticians have made great attempts to into a system of differential equations. develop statistical approaches for estimating ODE para- In a series of allometric studies by West et al. [32-34], meters by modeling the structure of measurement errors a power relationship that universally exists between [35-42]. We implemented Ramsay et al.’s [41] penalized parts and the whole can be explained by two fundamen- spline method for estimating constant dynamic para- tal design principles in biophysics and biochemistry; i.e., meters in our genetic mapping. The problem for sys- all organisms tend to maximize their metabolic capacity tems mapping with ODE models is different from those by increasing surface areas for energy and material considered in current literature. First, systems mapping Wu et al. BMC Systems Biology 2011, 5:84 Page 3 of 11 http://www.biomedcentral.com/1752-0509/5/84 is constructed within a mixture-based framework leaf (Figure 1A) and root biomass trajectories (Figure because QTL genotypes that define the DE models are 1C) delineate reasonably the decay of leaf and root bio- missing. Second, systems mapping incorporates genoty- mass at a late stage of development due to senescence. pic data which are categorical or binary. These two As expected, stem biomass growth does not experience characteristics determine the high complexity of our sta- such a decay (Figure 1B) although growth at the late tistical model and computational algorithm used for sys- stage tends to be stationary. tems mapping. By scanning the genetic linkage map composed of 950 molecular markers located in 25 linkage groups, we detected two significant QTLs, one named as biomass1 Results QTL detection that resides between markers GMKF167a and We develop a new model for QTL mapping by treating GMKF167b and the second as biomass2 that resides trait formation as a dynamic system and further incor- between markers sat 274 and BE801128 (Additional file porating the design principles of the biological system 1, Figure S1). Using the maximum likelihood estimates into a statistical mapping framework (see Methods). of the curve parameters in ODE (1) whose standard The new model, named systems mapping in light of its errors were obtained by the parameter bootstrap [43] systems feature of phenotypic description, was used to (Table 1), we drew the growth trajectories of leaf, stem map QTLs for biomass partitioning in a mapping popu- and root biomass for two different genotypes at each lation of soybeans composed of 184 recombinant inbred QTL (Figure 2). The genetic effects of the QTLs dis- lines (RIL) derived from two cultivars, Kefeng No. 1 and played different temporal patterns for three organs. The Nannong 1138-2. For an RIL population, there are two QTLs were expressed more rapidly with time for the homozygous genotypes, one composed of the Kefeng stem than for the leaves and roots. At biomass1,the No. 1 alleles and the other composed of the Nannong alleles from parent Nannong 1138-2 increase leaf and 1138-2 alleles. Figure 1 illustrates the growth trajectories stem biomass growth (Figure 2A and 2B), whereas the of leaf, stem and root biomass for individual RILs. By alleles from parent Kefeng No. 1 increase root biomass using the system of ODEs (1) to fit growth trajectories growth (Figure 2C). This could be interpreted as the of leaf, stem and root biomass over time, we obtained a Nannong 1138-2 allele favoring carbon allocation to the mean curve for each trait. It can be seen that growth shoots at the expense of the roots but the Kefeng No1 trajectories can be well modeled by three interconnect- allele favoring carbon allocation to the roots at the ing ODEs (1) derived from coordination theory [27] and expense of the shoots. Likewise, the biomass2 alleles allometric scaling [32-34]. The model-fitted curves of from Nannong 1138-2 favor carbon allocation to the A B C 2 4 6 8 2 4 6 8 2 4 6 8 Time Time Time Figure 1 Growth trajectories of leaf (A), stem (B) and root biomass (C) measured at multiple time points in a growing season of soybeans. Each grey line presents the growth trajectory of one of 184 RILs, whereas black lines are the mean growth trajectories of all RILs fitted using a system of ODEs (1). Biomass Wu et al. BMC Systems Biology 2011, 5:84 Page 4 of 11 http://www.biomedcentral.com/1752-0509/5/84 Table 1 The maximum likelihood point estimates (PEs) of ODE parameters and standard errors (SEs) of the estimates for the QTLs detected QTL Model Genotype Estimate a b g a b a b g L L L S S R R R 1M1 QQ PE 2.09 0.16 0.43 0.93 0.07 1.52 0.66 2.91 SE 0.02 1e-3 3e-3 0.01 2e-4 0.01 6e-3 0.02 qq PE 2.53 0.11 0.36 0.92 0.04 1.57 0.54 3.90 SE 0.01 4e-4 2e-3 0.01 1e-4 0.01 5e-3 0.01 M0 PE 2.30 0.13 0.39 0.92 0.05 1.56 0.60 3.37 SE 0.07 2e-3 0.01 0.02 4e-4 0.05 0.02 0.08 2M1 QQ PE 1.89 0.14 0.44 1.04 0.07 1.11 0.56 1.85 SE 0.04 1e-3 0.01 0.01 1e-4 0.01 5e-3 0.01 qq PE 2.55 0.10 0.31 0.98 0.04 1.11 0.51 2.18 SE 0.04 1e-3 0.01 0.01 4e-4 0.01 5e-3 0.02 M0 PE 2.25 0.12 0.37 1.03 0.05 1.12 0.55 2.06 SE 0.08 3e-3 0.02 0.02 6e-4 0.02 0.01 0.05 Model M1 assumes two different genotypes at a QTL (under the H ), whereas Model M0 assumes a single genotype (under the H ). 1 0 Note: QTL 1 is one between markers GMKF167a and GMKF167b on linkage group 12. QTL 2 is one between markers sat 274 and BE801128 on linkage group 24. The alleles of genotype QQ are derived from Kefeng No. 1 and those of genotype qq derived from Nannong 1138-2. leaves (Figure 2D) and those from Kefeng No. 1 favor Figure 4 shows the dynamic pattern of biomass parti- carbon allocation to the roots, but the alleles at this tioning to different organs. In general, the stem receives QTL inherited from parent Kefeng No. 1 favor carbon increasing allocation with time, whereas the partitioning allocation to the stem, which is different from the beha- to the leaves and roots decreases with time. Both QTLs vior of QTL biomass1. Note that leaf and root biomass detected, biomass1 and biomass2, control the degree of growth tend to decay at the late stage for almost all such time-dependent increase or decrease. For example, RILs. But the genotypes at the QTLs detected do not at biomass1, the Kefeng No. 1 genotype always exhibits reflect this trend (Figure 2), although they display much a larger degree of increasing biomass partitioning to the reduced rates of growth at the late stage. We explained stem but a larger degree of decreasing biomass parti- this to be due to some other QTLs that have not been tioning to the leaves and roots than the Nannong 1138- detected with the current linkage map. 2 genotype (Figure 4A vs. 4B). QTL biomass2 has a similar pattern of biomass partitioning for the stem and The functional relationships among leaf, stem and root biomass were determined by the QTLs detected (Fig- leaves, although it displays a stronger effect than does ure 3). For biomass1, two genotypes are not only dif- QTL biomass1.AtQTL biomass2,there is alarger ferent in whole-plant biomass trajectory, but also degree of decreasing biomass partitioning to the roots display pronounced discrepancies in biomass growth for the Nannong 1138-2 genotype than the Kefeng trajectories of individual organs (Figure 3A and 3B). No. 1 genotype (Figure 4C vs. 4D). This means that this QTL affects the dynamics of both plant size and biomass partitioning. The genotype Simulation composed of the Kefeng No. 1 alleles has a smaller By analyzing a real data set for soybean mapping, sys- slope of biomass growth, leading to smaller whole- tems mapping produces the identification of two signifi- plant biomass at late stages of development, than that cant QTLs that control the dynamic formation of composed of the Nannong 1138-2 alleles, but the for- whole-plant biomass through developmental regulation mer has larger root biomass over the entire period of of different organs, stem, leaves, and roots. To validate growth at the expense of the shoots than the latter. the new model, we performed simulation studies by For biomass2, two genotypes are similar in total plant mimicking the effects of QTL biomass2 detected from size during growth, but they have a marked distinction the example of QTL mapping in soybeans. The simu- in biomass partitioning (Figure 3C and 3D). It appears lated mapping population contains the same genotype that this QTL affects plant growth trajectories through data for 184 RILs. The phenotypic values of three traits, altering biomass partitioning rather than total amount the stem, leaf and root biomass, assumed to obey the of biomass. At this QTL, the genotype with the Kefeng system of ODE (1), were simulated at six different time No. 1 alleles has a dominant main stem and heavy points by summing time-dependent genotypic values at roots, whereas the genotype with the Nannong 1138-2 biomass2 calculated with curve parameters in Table 1 alleles carries dense leaves. and residual errors. Specifically, the phenotypic values of Wu et al. BMC Systems Biology 2011, 5:84 Page 5 of 11 http://www.biomedcentral.com/1752-0509/5/84 A: Leaf B: Stem C: Root Kefeng No. 1 Kefeng No. 1 Kefeng No. 1 Nannong 1138-2 Nannong 1138-2 Nannong 1138-2 2 4 6 8 2 4 6 8 2 4 6 8 D: Leaf E: Stem F: Root Kefeng No. 1 Kefeng No. 1 Kefeng No. 1 Nannong 1138-2 Nannong 1138-2 Nannong 1138-2 2 4 6 8 2 4 6 8 2 4 6 8 Time Time Time Figure 2 Growth trajectories of leaf (A, D), stem (B, E) and root biomass (C, F) for two different genotypes (presented by solid and broken black curves) at a QTL detected on linkage group 12 (upper panel) and 24 (lower panel), respectively. Two genotypes at a QTL are the homozygote for the alleles inherited from Kefeng No.1 (solid) and the homozygote for the allele from Nannong 1138-2 (broken). Curves in grey are growth trajectories of 184 RILs. the kth trait were simulated by adding white noise with of noise variance σ were set as the estimates from the variances σ to the ODE curves for the jth QTL geno- 2 2 k real data, which are σ =2.42, σ =1.72 and 1 2 type with the probability w , i.e., the conditional prob- j|i 2 σ =0.14 for the leaf, stem, and root biomass, respec- ability of the ith RIL that carries the jth QTL genotype, tively. Meanwhile, by assuming a modest heritability given the two markers genotypes of this RIL. The values Biomass Biomass Wu et al. BMC Systems Biology 2011, 5:84 Page 6 of 11 http://www.biomedcentral.com/1752-0509/5/84 Kefeng No.1 Alleles Nannong 1138-2 Alleles Whole Whole Leaf Leaf Stem Stem Root Root 2 4 6 8 2 4 6 8 C D Whole Whole Leaf Leaf Stem Stem Root Root 2 4 6 8 2 4 6 8 Time Time Figure 3 Growth trajectories of whole-plant (red), leaf (green), stem (blue) and root biomass (black) for two different genotypes at a QTL detected on linkage group 12 (upper panel) and 24 (lower panel). Two genotypes at a QTL are the homozygote for the alleles inherited from Kefeng No.1 (A, C) and the homozygote for alleles inherited from Nannong 1138-2 (B, D). (0.05) for each trait at a middle stage of growth, we re- simulated data is a mimicry of the real soybean data, the results suggest that the experimental design used for scaled σ values which were used to simulate a new soybean mapping is scientifically sound and can provide data set. convincing QTL detection. This actually is not surpris- Systems mapping, implemented with the parameter ing because phenotyping has low measurement errors. cascading method, estimates QTL genotype-specific In analyzing a simulated data set for the traits curve parameters in the ODE (1) from the simulated assumed to have a modest heritability (0.05), the esti- data. The simulation was repeated 100 times to calculate mates of the ODE parameters are reasonably accurate the means, biases, standard deviations, and root mean and precise, indicating the power of systems mapping to square errors, with results tabulated in Table 2. It was detect small QTLs involved in trait formation. We per- found that the model can provide reasonably accurate formed an additional simulation to investigate the and precise estimates of QTL genotype-specific ODE power and false positive rates of the model by changing parameters with a modest sample size (n = 184). The levels of noises. In general, the power of the model is biases of the estimates are negligible, compared with the high, reaching 0.80 even when the heritability of growth scale of the standard deviations. Given that this Biomass Biomass Wu et al. BMC Systems Biology 2011, 5:84 Page 7 of 11 http://www.biomedcentral.com/1752-0509/5/84 Kefeng No.1 Alleles Nannong 1138-2 Alleles Leaf Stem Root Leaf Stem Root 2 4 6 8 2 4 6 8 C D Leaf Stem Root Leaf Stem Root 2 4 6 8 2 4 6 8 Time Time Figure 4 Time-dependent percentages of leaf (green), stem (blue) and root biomass (thin black) for two different genotypes at a QTL detected on linkage group 12 (upper panel) and 24 (lower panel). Two genotypes at a QTL are the homozygote for the alleles inherited from Kefeng No.1 (A, D) and the homozygote for alleles from Nannong 1138-2 (B, C). curves is low (0.05). Basically, a QTL can be fully statistical model called functional mapping has been detected when the heritability is 0.10 or larger. In any developed by incorporating the mathematical aspects of case, the false positive rates are not beyond 0.10, mostly trait development [19-21]. Despite significant improve- being less than 0.05. ment over conventional static mapping, functional map- ping has still a major limitation in characterizing Discussion developmental pathways that cause a final phenotype Mapping the genetic architecture of complex traits has and unraveling the underlying genetic mechanisms for been a subject of interest in both theoretical and empiri- trait formation and progression. cal aspects of modern biology [3-15]. Original In this article, we have for the first time put forward a approaches for genetic mapping are based on single- new approach-systems mapping by treating a phenotypic point variation in a phenotypic trait, neglecting the trait as a dynamic system and incorporating the design dynamic change of the trait during development. To principles of a biological system into a statistical frame- capture the dynamic pattern of genetic control, a new work for genetic mapping. Various components that Biomass Percentage Biomass Percentage Wu et al. BMC Systems Biology 2011, 5:84 Page 8 of 11 http://www.biomedcentral.com/1752-0509/5/84 Table 2 Means of maximum likelihood estimates of curve reveal the genetic control mechanisms for several parameters from the ODE system (1) and their biases, mechanistically meaningful relationships. They are (1) standard deviations (STD) and square root mean square size-shape relationship - is a big plant due to a big stem errors (RMSE) from 100 simulation replicates with sparse leaves or a small stem with dense leaves? (2) Genotype Estimate a b g a b a b g L L L S S R R R structural-functional relationship - in a specific environ- QQ TRUE 2.55 0.10 0.31 0.98 0.04 1.11 0.51 2.18 ment does a plant tend to allocate more carbon to its MEAN 2.55 0.10 0.31 0.98 0.04 1.11 0.51 2.18 leaves for CO uptake or roots for water and nutrient BIAS*10 2.10 -0.07 0.72 0.27 0.01 3.49 -0.88 -2.89 uptake? (3) cause-effect relationship - are more roots STD*10 3.75 0.08 0.60 0.94 0.01 1.29 0.49 1.11 duetomoreleavesordomoreleavesproducemore RMSE*10 3.75 0.08 0.61 0.94 0.01 1.34 0.50 1.15 roots? and (4) pleiotropic relationship - different traits with a similar function tend to integrate into modularity qq TRUE 1.89 0.14 0.44 1.04 0.07 1.11 0.56 1.85 [45]. How do the same QTLs pleiotropically control this MEAN 1.89 0.14 0.44 1.04 0.07 1.11 0.56 1.85 modularity? A better understanding of these relation- BIAS*10 -4.23 0.11 0.55 -0.61 -0.05 1.62 0.92 -2.02 ships helps to gain more insights into the mechanistic STD*10 3.60 0.14 0.92 0.99 0.04 1.38 0.45 1.66 response of plant size and shape to developmental and RMSE*10 3.62 0.14 0.92 0.99 0.04 1.39 0.46 1.67 environmental signals and, also, provide guidance to select an ideotype of crop cultivars with optimal shape and structure suited to a particular environment [46]. constitute the system through developmental regulation The model described in this article is a simple frame- are studied and connected by a system of biologically work for systems mapping. It can be used as a start meaningful differential equations (DE). Thus, the genetic mapping of a complex phenotype become an issue of point to expand the concept of systems mapping to testing and estimating genotype-specific curve para- tackle more complicated biological problems. A pheno- meters at specific QTLs that define the emerging prop- type can be dissected to any number of components at erties and dynamic behavior of the DE system. Systems any level of organization, molecule, cell, tissue, or whole mapping identifies QTLs that control developmental organism, depending on the interest of researchers and interactions of traits, the temporal pattern of QTLs data availability. With more knowledge about phenotype expression during development, as well as the genetic formation and development, more components can be determinants that control developmental switches (on/ involved in a system that is specified by high-dimension off). differential equations. Sophisticated mathematical tech- Systems mapping was applied to map QTLs for niques are needed to obtain stable solutions of these equations. In addition, by integrating it with genome- dynamic trajectories of biomass from different organs, wide association studies, systems mapping will not only the stem, leaves, and roots, that interact and coordinate provide a clear view of how different components inter- to determine whole-plant biomass growth, in an experi- act and coordinate to form a phenotype, but also will be mental cross of soybeans between Kefeng No. 1 and capable of illustrating a comprehensive picture of the Nannong 1138-2. In general, the stem receives a propor- genetic architecture of complex phenotypes. There is tionally larger amount of biomass with development, also a good reason to integrate systems mapping with accompanying a proportional decrease of biomass to the network biology to explore how “omics” information leaves and roots. Specific QTLs, biomass1 and biomass2, contribute to the regulatory mechanisms of phenotype that control this allometric change with development formation [47]. In any case, systems mapping will open have been detected from systems mapping. The alleles a new avenue for understanding the genetic architecture at the two QTLs inherited from Kefeng No. 1 tend to of complex phenotypes from a perspective of mechanis- amplify this contrast in development-dependent biomass partitioning, as compared to those from Nannong 1138- tic pathways inside their formation. 2. One of the two QTLs, biomass2, was found in a simi- lar genomic region identified by traditional functional Conclusions mapping [44]. This consistency does not only simply The past two decades have seen a phenomenal increase verify our systems mapping, but also gains new insight in the number of tools for the genetic mapping of com- into biological functions of the detected QTLs. For plex traits. Although genetic mapping continues to be example, the two QTLs detected, biomass1 and bio- an interesting area in genetic research owing to the suc- mass2, trigger genetic effects on the interactions and cess of molecular and sequencing technologies in gener- coordination of different organs which cause the ating a flood of data, a conceptual breakthrough in this dynamic variation of biomass growth. area remains elusive. In this article, we present a bot- Through various tests for ODE parameters individu- tom-top model for mapping and studying the genetic ally or in a combination, our systems mapping can architecture of complex traits. Different from existing Wu et al. BMC Systems Biology 2011, 5:84 Page 9 of 11 http://www.biomedcentral.com/1752-0509/5/84 mapping models, we use a systems approach to identify errors for such a large-scale field trial. Phenotyping pre- specific genes or quantitative trait loci that govern the cision was estimated to be above 95%. developmental interactions of various components com- Unlike a traditional mappingproject,ourgoal is to prising the phenotype. The map of developmental inter- map QTLs that control the dynamic process of how dif- actions among different components is constructed by a ferent organs, the stem, leaves, and roots, interact and system of differential equations. Thus, by estimating and coordinate to determine whole-plant biomass. The inter- testing mathematical parameters that specify the system, actions and coordination of different organs for a plant are understood using design principles described by the we are able to predict or alter the physiological status of ODE system (1). a phenotype based on the underlying genetic control mechanisms. We have tested and validated our model by analyzing a real data set for genetic mapping of bio- Likelihood mass growth in soybeans. The detection of QTLs by the Let denote the vector of y = y (t ) , ··· , y t ki i1 ki im ki i new model provides biologically meaningful interpreta- phenotypic values for trait k (k = 1 for leaf biomass (L), tions of QTL effects on trait formation and dynamics. 2 for stem biomass (S), and 3 for root biomass (R)) mea- The new model can be readily used to study the genetic sured on progeny i at time points .Note t , ··· , t i1 im basis of phenotypes in any other organism. i that the number of time points measured may be pro- geny-specific, expressed as m for progeny i.Assuming Methods that multiple QTLs (segregating with J genotypes), each Mapping Population bracketed by two flanking markers M, affects these Our model derivation is based on a mapping population three traits, we construct a mixture model-based likeli- comprising of n recombinant inbred lines (RILs), hood as initiated with two inbred lines. By continuous selfing or n J inbreeding, RILs after the F generation are considered L(z, M)= [ω f (z ;  , )] (2) homozygous because the fixation at any locus is given j|i j i j 7-1 i=1 j=1 by f = 1 - 0.5 ≈ 1. In practical terms, all plants from a single RIL are genetically identical, and can be used for where y =(y , y , y ) is a joint vector of phenotypic 1 2 3 replicated experiments under different environments. In values for the three traits, with z =(y , y , y )present- i 1i 2i 3i addition, each RIL represents a unique combination of ing the z-vector for progeny i; ω is the conditional j|i alleles from the parental genotypes where there are two probability of QTL genotype j (j = 1,..., J) given the mar- homozygous genotypes at each marker locus, each cor- ker genotype of progeny i, which can be expressed as a responding to a parental allele. The mapping population function of the recombination fractions between the is genotyped at molecular markers to construct a linkage QTL and markers [50], and f (z ; Θ , Ψ)isanMVN of j i j map covering the entire genome. The recombination leaf, stem and root biomass for progeny i which carries fraction between two markers is converted to the QTL genotype j, with mean vectors genetic distance in centiMorgan (cM) through a map function, such as the Haldane or Kosambi map function. μ =(μ , μ , μ ), j =1, ... , J j 1j 2j 3j The map constructed is used to locate QTLs that con- trol a quantitative trait of interest. specified by Θ , and covariance matrix specified by Ψ. We obtained a sample of 184 RILs derived from two If a system of differential equations (1) is used to jointly cultivars, Kefeng No. 1 and Nannong 1138-2, for map- model QTL genotype-specific means vectors for the ping agronomic traits. These RILs were genotyped for three traits, then we have Θ=(a , b , l , a , b , a , j Lj Lj jL Sj Sj Rj 950 molecular markers locating in 25 linkage groups b , l ) for genotype j. Rj Rj [48,49]. The plants were grown in a simple lattice design with two replicates in a plot at Jiangpu Soybean Experi- Estimation ment Station, Nanjing Agricultural University, China. Unlike a traditional mixture model for QTL mapping, Ten plants in the second row of a plot were randomly we will model the genotypic values of each QTL geno- type in likelihood (2) characterized by a group of non- selected for measuring leaf, stem and root biomass at linear ODEs. While an analytical solution is not each time in the whole growing season. After 20 days of available, we will implement numerical approaches to seedling emergence, dry weights separately for the solve these ODEs. Let μ (t) denote the genotypic value leaves, stem and roots were measured once every 5 to kj of the kth trait at time t for a QTL genotype j.Thus, 10 days until most plants stopped growth. A total of 6 the dynamic system of the traits and their interactions, to 8 measurements were taken for each of the RILs stu- died. Great efforts were made to control measurement regulated by QTL genotypes, can be modeled by a Wu et al. BMC Systems Biology 2011, 5:84 Page 10 of 11 http://www.biomedcentral.com/1752-0509/5/84 system of ODE (1), Additional material dμ (t) kj Additional file 1: Figure S1. The profiles of the log-likelihood ratios (LR) (3) = g (μ (t),  ), k =1, ··· ,3, k j between the full model (there is a QTL) and reduced model (there is no dt QTL) for soybean height growth trajectories throughout the soybean genome composed of 25 linkage groups. where μ (t)= (μ (t), ..., μ (t)) ,and Θ is a vector of j 1j 3j j ODE parameters associated with QTL genotype j.For J possible genotypes in the mapping population, we have T T Acknowledgements . The question now is how to Θ esti- =  , ... , 1 J This work is supported by the Changjiang Scholars Award, “One-thousand mate from noisy measurements. The functional mean Person” Award, NSF/IOS-0923975, a discovery grant of the Natural Sciences and Engineering Research Council of Canada (NSERC) (J. Cao), the National μ (t) may be represented as a linear combination of kj Key Basic Research Program of China (2009CB1184，2010CB1259， basis functions: 2011CB1093), The National High Technology R&D Program of China (2009AA1011), and the MOE 111 Project (B08025). Part of this work was carried when RW and JC were invited Research Fellows at the Statistical and μ (t)= c φ (t)= c φ (t) (4) kj kjr kjr Applied Mathematical Sciences Institute (SAMSI), sponsored by Duke kj kj University, University of North Carolina at Chapel Hill, and North Carolina r=1 State University. Thanks are due to Prof. Shouyi Chen and Prof. Deyue Yu for kind permission to use the jointly developed NJRIKY genetic linkage map. where j (t) = (j (t), ..., j (t)) is a vector of basis kj kj1 kjR functions with R orders and c =(c ,..., c ) is a vec- kj kj1 kjR Author details 3;J Center for Computational Biology, National Engineering Laboratory for Tree tor of basis coefficients. Define c = c as a kj k=1;j=1 Breeding, Key Laboratory of Genetics and Breeding in Forest Trees and Ornamental Plants, Beijing Forestry University, Beijing 100083, China. length (R ×3× J) vector of basis coefficients. The cubic Department of Statistics & Actuarial Science, Simon Fraser University, B-splines are often chosen as basis functions, since any Burnaby, B.C. Canada V5A 1S6. Department of Agronomy, Henan Institute of B-spline basis function is only positive over a short sub- Science and Technology, Xinxiang 453003, China. Center for Statistical Genetics, Pennsylvania State University, Hershey, PA 17033, USA. National interval and zero elsewhere. This is called the compact Center for Soybean Improvement, National Key Laboratory for Crop Genetics support property, and is essential for efficient computa- and Germplasm Enhancement, Soybean Research Institute, Nanjing tion. The flexibility of the B-spline basis functions Agricultural University, Nanjing 210095, China. Department of Horticultural Sciences, University of Florida, Gainesville, FL 32611, USA. depend on the number and location of knots we choose. It is an infinite-dimension optimization problem to Authors’ contributions choose the optimal number of knots and their locations. RW conceived of the idea of systems mapping, coordinated the whole study, and wrote the manuscript. JC derived the model, computed the real A popular approach to avoid this dilemma is choosing a data, run simulation studies, and participated in the writing of the Methods saturated number of knots and using a roughness pen- section. ZH conducted the soybean experiment and collected phenotypic alty to control the smoothness of the fitted curve and data. ZW packed the model into a computer package SysMap. JG directed the experimental design and data collection of soybeans. EV oversaw the avoid over-fitting [40]. project and highlighted the biological relevance a computational model We estimate the basis coefficient c and ODE para- must possess. All authors read and approved the final manuscript. meter Θ based on a two-nested level of optimization. In Received: 7 February 2011 Accepted: 27 May 2011 the inner level of optimization, c is estimated by opti- Published: 27 May 2011 mizing a criterion U(c|Θ), given any value of Θ.There- fore, the estimate maybeviewedasa function of Θ, cˆ References which is denoted as ˆ . Since no analytic formula 1. Lee SH, van der Werf JHJ, Hayes BJ, Goddard ME, Visscher PM: Predicting c () unobserved phenotypes for complex traits from whole-genome SNP for is available, is an implicit function. In the cˆ () cˆ data. PLoS Genet 2008, 4(10):e1000231.. outer level of optimization, Θ is estimated by optimizing 2. Lynch M, Walsh B: Genetics and Analysis of Quantitative Traits Sinauer Associates, Sunderland, MA; 1998. a criterion H , cˆ () . The parameter is removed in cˆ 3. Lander ES, Botstein D: Mapping Mendelian factors underlying quantitative traits using RFLP linkage maps. Genetics 1989, 121:185-199. the parameter space in the outer level by treating it as 4. Zeng ZB: Precision mapping of quantitative trait loci. Genetics 1994, an implicit function of Θ. Although cˆ () does not 136:1457-1468. have an analytic formula, the outer level of optimization 5. Xu S, Atchley WR: A random model approach to interval mapping of quantitative genes. Genetics 1995, 141:1189-1197. only requires to calculate the derivative dcˆ/d ,which 6. Wu RL, Ma CX, Casella G: Joint linkage and linkage disequilibrium can be obtained by using the implicit function theorem. mapping of quantitative trait loci in natural populations. Genetics 2002, 160:779-792. The above optimization procedure is called the para- 7. Yi N, Xu S: Bayesian lasso for quantitative trait loci mapping. Genetics meter cascading method. Note that when the two cri- 2008, 179:1045-1055. teria U(c|Θ)and H , cˆ () are the same, the 8. Zou F, Nie L, Wright FA, Sen PK: A robust QTL mapping procedure. J Stat Plann Infer 2009, 139:978-989. parameter cascading method is equivalent to the profil- 9. Ehrenreich IM, Torabi N, Jia Y, Kent J, Martis S, Shapiro JA, Gresham D, ing method. Caudy AA, Kruglyak L: Dissection of genetically complex traits with extremely large pools of yeast segregants. Nature 2010, 464:1039-1042. Wu et al. BMC Systems Biology 2011, 5:84 Page 11 of 11 http://www.biomedcentral.com/1752-0509/5/84 10. Steinmetz LM, Sinha H, Richards DR, Spiegelman JI, Oefner PJ, McCusker JH, 37. Huang JZ, Liu N, Pourahmadi M, Liu L: Covariance selection and Davis RW: Dissecting the complex architecture of a quantitative trait estimation via penalized normal likelihood. Biometrika 2006, 93:85-98. locus in yeast. Nature 2002, 416:326-330. 38. Li L, Brown MB, Lee KH, Gupta S: Estimation and inference for a spline- 11. Singer JB, Hill AE, Burrage LC, Olszens KR, Song J, Justice M, O’Brien WE, enhanced population pharmacokinetic model. Biometrics 2002, Conti DV, Witte JS, Lander ES, Nadeau JH: Genetic dissection of complex 58:601-611. traits with chromosome substitution strains of mice. Science 2004, 39. Ramsay JO: Principal differential analysis: Data reduction by differential 304:445-448. operators. J Roy Stat Soc Ser B 1996, 58:495-508. 12. Altshuler D, Daly MJ, Lander ES: Genetic mapping in human disease. 40. Ramsay JO, Silverman BW: Functional Data Analysis. 2 edition. Springer, New Science 2008, 322:881-888. York; 2005. 13. Shao H, Burrage LC, Sinasac DS, Hill AE, Ernest SR, O’Brien W, Courtland HW, 41. Ramsay JO, Hooker G, Campbell D, Cao JG: Parameter estimation for Jepsen KJ, Kirby A, Kulbokas EJ, Daly MJ, Broman KW, Lander ES, Nadeau JH: differential equations: a generalized smoothing approach (with Genetic architecture of complex traits: Large phenotypic effects and discussion). J Roy Stat Soc Ser B 2007, 69:741-796. pervasive epistasis. Proc Natl Acad Sci USA 2008, 105:19910-19914. 42. Liang H, Wu HL: Parameter estimation for differential equation models 14. Mackay TFC, Stone EA, Ayroles JF: The genetics of quantitative traits: using a frame-work of measurement error in regression model. JAm challenges and prospects. Nat Rev Genet 2009, 10:565-577. Stat Assoc 2008, 103:1570-1583. 15. Balasubramanian S, Schwartz C, Singh A, Warthmann N, Kim MC, Maloof JN, 43. Efron B, Tibshirani RJ: An Introduction to the Bootstrap Chapman & Hall; Loudet O, Trainer GT, Dabi T, Borevitz JO, Chory J, Weigel D: QTL mapping 1993. in new Arabidopsis thaliana advanced intercross-recombinant inbred 44. Li Q, Huang Z, Xu M, Wang C, Gai J, Huang Y, Pang X, Wu RL: Functional lines. PLoS ONE 2009, 4(2):e4318. mapping of genotype-environment interactions for soybean growth by 16. Frary A, Nesbitt TC, Grandillo S, Knaap E, Cong B, Liu J, Meller J, Elber R, a semiparametric approach. Plant Methods 2010, 6:13. Alpert KB, Tanksley SD: fw2. 2: a quantitative trait locus key to the 45. McCarthy MC, Enquist BJ: Consistency between an allometric approach evolution of tomato fruit size. Science 2000, 289:85-88. and optimal partitioning theory in global patterns of plant biomass 17. Li C, Zhou A, Sang T: Rice domestication by reducing shattering. Science allocation. Funct Ecol 2007, 21:713-720. 2006, 311:1936-1939. 46. Jiao Y, Wang Y, Xue D, Wang J, Yan M, Liu G, Dong G, Zeng D, Lu Z, Zhu X, 18. Huang X, Qian Q, Liu Z, Sun H, He S, Luo D, Xia G, Chu C, Li J, Fu X: Qian Q, Li J: Regulation of OsSPL14 by OsmiR156 defines ideal plant Natural variation at the DEP1 locus enhances grain yield in rice. Nat architecture in rice. Nat Genet 2010, 42:541-544. Genet 2009, 41:494-497. 47. Wu S, Yap JS, Li Y, Li Q, Fu GF, Li JH, Das K, Berg A, Zeng YR, Wu RL: 19. Ma CX, Casella G, Wu RL: Functional mapping of quantitative trait loci Network models for dissecting plant development by functional underlying the character process: a theoretical framework. Genetics 2002, mapping. Curr Bioinform 2009, 4:183-187. 161:1751-1762. 48. Li H, Huang Z, Gai J, Wu S, Zeng Y, Wu RL: A conceptual framework for 20. Wu RL, Lin M: Functional mapping - how to map and study the genetic mapping quantitative trait loci regulating ontogenetic allometry. PLoS architecture of dynamic complex traits. Nat Rev Genet 2006, 7:229-237. ONE 2007, 2(11):e1245. 21. Li Y, Wu RL: Functional mapping of growth and development. Biol Rev 49. Zhang W-K, Wang Y-J, Luo G-Z, Zhang J-S, He C-Y, Wu XL, Gai JY, Chen SY: 2010, 85:207-216. QTL mapping of ten agronomic traits on the soybean (Glycine max L. 22. Kitano H: Systems biology: a brief overview. Science 2002, 295:1662-1664. Merr.) genetic map and their association with EST markers. Theor Appl 23. Jansen RC: Studying complex biological systems using multifactorial Genet 2004, 108:1131-1139. perturbation. Nat Rev Genet 2003, 4:145-151. 50. Wu RL, Ma CX, Casella G: Statistical Genetics of Quantitative Traits: Linkage, 24. Csete ME, Doyle JC: Reverse engineering of biological complexity. Science Maps, and QTL Springer-Verlag, New York; 2007. 2002, 295:1664-1669. doi:10.1186/1752-0509-5-84 25. Cannell MGR, Dewar RC: Carbon allocation in trees: a review of concepts Cite this article as: Wu et al.: Systems mapping: how to improve the for modelling. Ad Ecol Res 1994, 25:59-104. genetic mapping of complex traits through design principles of 26. Luo Y, Field CB, Mooney HA: Predicting responses of photosynthesis and biological systems. BMC Systems Biology 2011 5:84. root fraction to elevated CO Interactions among carbon, nitrogen, and 2: growth. Plant Cell Environ 1994, 17:1195-1204. 27. Chen J, Reynolds J: A coordination model of carbon allocation in relation to water supply. Ann Bot 1997, 80:45-55. 28. Weiner J: Allocation, plasticity and allometry in plants. Perspect Plant Ecol Evol Syst 2004, 6:207-215. 29. Hermans C, Hammond JP, White PJ, Verbruggen N: How do plants respond to nutrient shortage by biomass allocation? Trends Plant Sci 2006, 11:610-617. 30. Marcelis LFM, Heuvelink E: Concepts of modelling carbon allocation among plant organs. In Functional-Structural Plant Modelling in Crop Production. Edited by: Vos J, Marcelis LFM, de Visser PHB, Struik PC and Evers JB. Springer. Printed in the Netherlands; 2007:103-111. 31. Genard M, Dauzat J, Franck N, Lescourret F, Moitrier N, Vaast P, Vercambre G: Carbon allocation in fruit trees: From theory to modelling. Trees - Structure and Function 2008, 22:269-282. 32. West GB, Brown JH, Enquist BJ: A general model for the origin of Submit your next manuscript to BioMed Central allometric scaling laws in biology. Science 1997, 276:122-126. and take full advantage of: 33. West GB, Brown JH, Enquist BJ: The fourth dimension of life: Fractal geometry and allometric scaling of organisms. Science 1999, • Convenient online submission 284:1677-1679. 34. West GB, Brown JH, Enquist BJ: A general model for ontogenetic growth. • Thorough peer review Nature 2001, 413:628-631. • No space constraints or color figure charges 35. Putter H, Heisterkamp SH, Lange JMA, De Wolf F: A Bayesian approach to • Immediate publication on acceptance parameter estimation in HIV dynamical models. Stat Med 2002, 21:2199-2214. • Inclusion in PubMed, CAS, Scopus and Google Scholar 36. Huang Y, Wu H: A Bayesian approach for estimating antiviral efficacy in • Research which is freely available for redistribution HIV dynamic models. J Appl Stat 2006, 33:155-174. Submit your manuscript at www.biomedcentral.com/submit http://www.deepdyve.com/assets/images/DeepDyve-Logo-lg.png BMC Systems Biology Springer Journals

Systems mapping: how to improve the genetic mapping of complex traits through design principles of biological systems

Download PDF
11 pages

Loading next page...
 
/lp/springer-journals/systems-mapping-how-to-improve-the-genetic-mapping-of-complex-traits-upa3g0rLm7

References (113)

Publisher
Springer Journals
Copyright
Copyright © 2011 by Wu et al; licensee BioMed Central Ltd.
Subject
Life Sciences; Bioinformatics; Systems Biology; Simulation and Modeling; Computational Biology/Bioinformatics; Physiological, Cellular and Medical Topics; Algorithms
eISSN
1752-0509
DOI
10.1186/1752-0509-5-84
pmid
21615967
Publisher site
See Article on Publisher Site

Abstract

Background: Every phenotypic trait can be viewed as a “system” in which a group of interconnected components function synergistically to yield a unified whole. Once a system’s components and their interactions have been delineated according to biological principles, we can manipulate and engineer functionally relevant components to produce a desirable system phenotype. Results: We describe a conceptual framework for mapping quantitative trait loci (QTLs) that control complex traits by treating trait formation as a dynamic system. This framework, called systems mapping, incorporates a system of differential equations that quantifies how alterations of different components lead to the global change of trait development and function through genes, and provides a quantitative and testable platform for assessing the interplay between gene action and development. We applied systems mapping to analyze biomass growth data in a mapping population of soybeans and identified specific loci that are responsible for the dynamics of biomass partitioning to leaves, stem, and roots. Conclusions: We show that systems mapping implemented by design principles of biological systems is quite versatile for deciphering the genetic machineries for size-shape, structural-functional, sink-source and pleiotropic relationships underlying plant physiology and development. Systems mapping should enable geneticists to shed light on the genetic complexity of any biological system in plants and other organisms and predict its physiological and pathological states. Background success of this prediction depends on how well we can Predicting the phenotype from the genotype of complex map the underlying QTLs and characterize complex organisms is one of the most important and challenging interactions of these QTLs with each other and with questions we face in modern biology and medicine [1]. environmental factors. Powerful statistical models have Genetic mapping, dissecting a phenotypic trait to its been developed in the past two decades to detect QTLs underlying quantitative trait loci (QTLs) through the and study their biological function in a diverse array of use of molecular markers, has proven powerful for phenotypic traits [3-9]. Worldwide, a substantial effort establishing genotype-phenotype relationships and pre- has been made resulting in the collection of a large dicting phenotypes of individual organisms based on amount of data aimed at the identification of QTLs [10-15]. Unfortunately, despite hundreds of thousands of their QTL genotypes responsible for the trait [2]. The QTLs detected in a diversity of organisms, only a few of them have been isolated by positional cloning (see * Correspondence: [email protected]; [email protected] [16-18]), leaving it unsolved how to construct a geno- † Contributed equally Center for Computational Biology, National Engineering Laboratory for Tree type-phenotype relationship map using genetic mapping. Breeding, Key Laboratory of Genetics and Breeding in Forest Trees and The most likely reason for this result may arise from a Ornamental Plants, Beijing Forestry University, Beijing 100083, China possibility that the QTLs detected by stringent statistical National Center for Soybean Improvement, National Key Laboratory for Crop Genetics and Germplasm Enhancement, Soybean Research Institute, tests are not biologically relevant. Existing strategies for Nanjing Agricultural University, Nanjing 210095, China QTL mapping were built on testing for a direct Full list of author information is available at the end of the article © 2011 Wu et al; licensee BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited. Wu et al. BMC Systems Biology 2011, 5:84 Page 2 of 11 http://www.biomedcentral.com/1752-0509/5/84 association between genotype and end-point phenotype. production as well as internal efficiency through redu- Although such strategies are simple and have been cing distances and the time to transport water, nutri- widely accepted, they neglect the biological processes ents, and carbon. The integration of this optimization involved in trait development [14]. To attempt to fill theory expressed in terms of allometric scaling with the this gap, a statistical model, called functional mapping coordination theory leads to a tripled group of ordinary [19-21], has been developed to study the interplay differential equations (ODEs) to specify the coordination of the leaf, stem, and root biomass for a plant: between genetics and the developmental process of a phenotypic trait by integrating mathematical models and dM computational algorithms. If a trait is understood as a = α W − γ M L L L dt “system” that is composed of many underlying biological dM components [22-24], we should be in a better position β (1) = α W dt to comprehend the process and behavior of trait forma- dM tion based on interactive relationships among different R = α W − γ M R R R components. Through mapping and using those QTLs dt that govern design principles of a biological system, a new trait that is able to maximize resource-use effi- where M , M ,and M are the biomasses of the leaves L S R ciency can be generated and engineered. (L), the stem (S), and the roots (R), respectively, with As one important strategy for plants to respond to whole-plant biomass W = M + M + M ; a and b are L S R variation in the availability of resources in their environ- the constant and exponent power of an organ biomass ment, biomass allocation has been extensively used to scaling as whole-plant biomass [32,33]; and g is the rate study the relationship between structure and function in of eliminating ageing leaves and roots. The complex modern ecology [25-28]. The concept of biomass alloca- interactions between different parts of a plant that tion has now been increasingly integrated with plant underlie design principles of plant biomass growth can management and breeding, aimed to direct a maximum be modeled and studied by estimating and testing the amount of biomass to the target of harvest (leaves, stem, ODE parameters (a , b , l , a , b , a , b , l ). For L L L S S R R R roots, or fruits) [29-31]. If the whole-plant biomass is example, plants are equipped with a capacity to optimize considered as a target trait, we need to understand how their fitness under low nutrient availability by shifting different organs of a plant coordinate and interact to the partitioning of carbohydrates to processes associated optimize the capture of nutrients, light, water, and car- with nutrient uptake at a cost of carbon acquisition bon dioxide in a manner that maximizes plant growth [29]. These parameters can be used to quantify and pre- rate through a specific developmental program because dict such regulation between different plant parts in plant biomass growth is not simply the addition of bio- response to environmental and developmental changes. mass to various organs. Many theories and models have In this article, we put forward a conceptual framework been proposed to predict the pattern of biomass parti- to incorporate the design principles of trait formation tioning in a response to changing environment. Chen and development into a statistical framework for QTL and Reynolds [27] used coordination theory to model mapping. Complementary to our previous functional the dynamic allocation of carbon to different organs mapping [19-21], we name this new mapping framework during growth in relation to carbon and water/nitrogen “systems mapping” in light of its systems dissection and supply by a group of differential equations. Compared modeling of phenotypic formation. A group of ODEs to the conventional optimization model in the context like (1) or other types of differential equations is used to of maximizing the relative growth rate of a plant, the quantify the phenotypic system. Much work in solving coordination model does not require an unrealistic ODEs has focused on the simulation and analysis of the capacity the plant possesses for knowing beforehand the behavior of state variables for a dynamic system, but the environmental conditions it will experience during the estimation of ODE parameters that define the system growth period. Here, we integrate the coordination and based on the measurement of state variables at multiple optimization model to study the pattern of biomass par- time points is relatively a new area. Yet, in the recent titioning by incorporating the allometric scaling theory years, many statisticians have made great attempts to into a system of differential equations. develop statistical approaches for estimating ODE para- In a series of allometric studies by West et al. [32-34], meters by modeling the structure of measurement errors a power relationship that universally exists between [35-42]. We implemented Ramsay et al.’s [41] penalized parts and the whole can be explained by two fundamen- spline method for estimating constant dynamic para- tal design principles in biophysics and biochemistry; i.e., meters in our genetic mapping. The problem for sys- all organisms tend to maximize their metabolic capacity tems mapping with ODE models is different from those by increasing surface areas for energy and material considered in current literature. First, systems mapping Wu et al. BMC Systems Biology 2011, 5:84 Page 3 of 11 http://www.biomedcentral.com/1752-0509/5/84 is constructed within a mixture-based framework leaf (Figure 1A) and root biomass trajectories (Figure because QTL genotypes that define the DE models are 1C) delineate reasonably the decay of leaf and root bio- missing. Second, systems mapping incorporates genoty- mass at a late stage of development due to senescence. pic data which are categorical or binary. These two As expected, stem biomass growth does not experience characteristics determine the high complexity of our sta- such a decay (Figure 1B) although growth at the late tistical model and computational algorithm used for sys- stage tends to be stationary. tems mapping. By scanning the genetic linkage map composed of 950 molecular markers located in 25 linkage groups, we detected two significant QTLs, one named as biomass1 Results QTL detection that resides between markers GMKF167a and We develop a new model for QTL mapping by treating GMKF167b and the second as biomass2 that resides trait formation as a dynamic system and further incor- between markers sat 274 and BE801128 (Additional file porating the design principles of the biological system 1, Figure S1). Using the maximum likelihood estimates into a statistical mapping framework (see Methods). of the curve parameters in ODE (1) whose standard The new model, named systems mapping in light of its errors were obtained by the parameter bootstrap [43] systems feature of phenotypic description, was used to (Table 1), we drew the growth trajectories of leaf, stem map QTLs for biomass partitioning in a mapping popu- and root biomass for two different genotypes at each lation of soybeans composed of 184 recombinant inbred QTL (Figure 2). The genetic effects of the QTLs dis- lines (RIL) derived from two cultivars, Kefeng No. 1 and played different temporal patterns for three organs. The Nannong 1138-2. For an RIL population, there are two QTLs were expressed more rapidly with time for the homozygous genotypes, one composed of the Kefeng stem than for the leaves and roots. At biomass1,the No. 1 alleles and the other composed of the Nannong alleles from parent Nannong 1138-2 increase leaf and 1138-2 alleles. Figure 1 illustrates the growth trajectories stem biomass growth (Figure 2A and 2B), whereas the of leaf, stem and root biomass for individual RILs. By alleles from parent Kefeng No. 1 increase root biomass using the system of ODEs (1) to fit growth trajectories growth (Figure 2C). This could be interpreted as the of leaf, stem and root biomass over time, we obtained a Nannong 1138-2 allele favoring carbon allocation to the mean curve for each trait. It can be seen that growth shoots at the expense of the roots but the Kefeng No1 trajectories can be well modeled by three interconnect- allele favoring carbon allocation to the roots at the ing ODEs (1) derived from coordination theory [27] and expense of the shoots. Likewise, the biomass2 alleles allometric scaling [32-34]. The model-fitted curves of from Nannong 1138-2 favor carbon allocation to the A B C 2 4 6 8 2 4 6 8 2 4 6 8 Time Time Time Figure 1 Growth trajectories of leaf (A), stem (B) and root biomass (C) measured at multiple time points in a growing season of soybeans. Each grey line presents the growth trajectory of one of 184 RILs, whereas black lines are the mean growth trajectories of all RILs fitted using a system of ODEs (1). Biomass Wu et al. BMC Systems Biology 2011, 5:84 Page 4 of 11 http://www.biomedcentral.com/1752-0509/5/84 Table 1 The maximum likelihood point estimates (PEs) of ODE parameters and standard errors (SEs) of the estimates for the QTLs detected QTL Model Genotype Estimate a b g a b a b g L L L S S R R R 1M1 QQ PE 2.09 0.16 0.43 0.93 0.07 1.52 0.66 2.91 SE 0.02 1e-3 3e-3 0.01 2e-4 0.01 6e-3 0.02 qq PE 2.53 0.11 0.36 0.92 0.04 1.57 0.54 3.90 SE 0.01 4e-4 2e-3 0.01 1e-4 0.01 5e-3 0.01 M0 PE 2.30 0.13 0.39 0.92 0.05 1.56 0.60 3.37 SE 0.07 2e-3 0.01 0.02 4e-4 0.05 0.02 0.08 2M1 QQ PE 1.89 0.14 0.44 1.04 0.07 1.11 0.56 1.85 SE 0.04 1e-3 0.01 0.01 1e-4 0.01 5e-3 0.01 qq PE 2.55 0.10 0.31 0.98 0.04 1.11 0.51 2.18 SE 0.04 1e-3 0.01 0.01 4e-4 0.01 5e-3 0.02 M0 PE 2.25 0.12 0.37 1.03 0.05 1.12 0.55 2.06 SE 0.08 3e-3 0.02 0.02 6e-4 0.02 0.01 0.05 Model M1 assumes two different genotypes at a QTL (under the H ), whereas Model M0 assumes a single genotype (under the H ). 1 0 Note: QTL 1 is one between markers GMKF167a and GMKF167b on linkage group 12. QTL 2 is one between markers sat 274 and BE801128 on linkage group 24. The alleles of genotype QQ are derived from Kefeng No. 1 and those of genotype qq derived from Nannong 1138-2. leaves (Figure 2D) and those from Kefeng No. 1 favor Figure 4 shows the dynamic pattern of biomass parti- carbon allocation to the roots, but the alleles at this tioning to different organs. In general, the stem receives QTL inherited from parent Kefeng No. 1 favor carbon increasing allocation with time, whereas the partitioning allocation to the stem, which is different from the beha- to the leaves and roots decreases with time. Both QTLs vior of QTL biomass1. Note that leaf and root biomass detected, biomass1 and biomass2, control the degree of growth tend to decay at the late stage for almost all such time-dependent increase or decrease. For example, RILs. But the genotypes at the QTLs detected do not at biomass1, the Kefeng No. 1 genotype always exhibits reflect this trend (Figure 2), although they display much a larger degree of increasing biomass partitioning to the reduced rates of growth at the late stage. We explained stem but a larger degree of decreasing biomass parti- this to be due to some other QTLs that have not been tioning to the leaves and roots than the Nannong 1138- detected with the current linkage map. 2 genotype (Figure 4A vs. 4B). QTL biomass2 has a similar pattern of biomass partitioning for the stem and The functional relationships among leaf, stem and root biomass were determined by the QTLs detected (Fig- leaves, although it displays a stronger effect than does ure 3). For biomass1, two genotypes are not only dif- QTL biomass1.AtQTL biomass2,there is alarger ferent in whole-plant biomass trajectory, but also degree of decreasing biomass partitioning to the roots display pronounced discrepancies in biomass growth for the Nannong 1138-2 genotype than the Kefeng trajectories of individual organs (Figure 3A and 3B). No. 1 genotype (Figure 4C vs. 4D). This means that this QTL affects the dynamics of both plant size and biomass partitioning. The genotype Simulation composed of the Kefeng No. 1 alleles has a smaller By analyzing a real data set for soybean mapping, sys- slope of biomass growth, leading to smaller whole- tems mapping produces the identification of two signifi- plant biomass at late stages of development, than that cant QTLs that control the dynamic formation of composed of the Nannong 1138-2 alleles, but the for- whole-plant biomass through developmental regulation mer has larger root biomass over the entire period of of different organs, stem, leaves, and roots. To validate growth at the expense of the shoots than the latter. the new model, we performed simulation studies by For biomass2, two genotypes are similar in total plant mimicking the effects of QTL biomass2 detected from size during growth, but they have a marked distinction the example of QTL mapping in soybeans. The simu- in biomass partitioning (Figure 3C and 3D). It appears lated mapping population contains the same genotype that this QTL affects plant growth trajectories through data for 184 RILs. The phenotypic values of three traits, altering biomass partitioning rather than total amount the stem, leaf and root biomass, assumed to obey the of biomass. At this QTL, the genotype with the Kefeng system of ODE (1), were simulated at six different time No. 1 alleles has a dominant main stem and heavy points by summing time-dependent genotypic values at roots, whereas the genotype with the Nannong 1138-2 biomass2 calculated with curve parameters in Table 1 alleles carries dense leaves. and residual errors. Specifically, the phenotypic values of Wu et al. BMC Systems Biology 2011, 5:84 Page 5 of 11 http://www.biomedcentral.com/1752-0509/5/84 A: Leaf B: Stem C: Root Kefeng No. 1 Kefeng No. 1 Kefeng No. 1 Nannong 1138-2 Nannong 1138-2 Nannong 1138-2 2 4 6 8 2 4 6 8 2 4 6 8 D: Leaf E: Stem F: Root Kefeng No. 1 Kefeng No. 1 Kefeng No. 1 Nannong 1138-2 Nannong 1138-2 Nannong 1138-2 2 4 6 8 2 4 6 8 2 4 6 8 Time Time Time Figure 2 Growth trajectories of leaf (A, D), stem (B, E) and root biomass (C, F) for two different genotypes (presented by solid and broken black curves) at a QTL detected on linkage group 12 (upper panel) and 24 (lower panel), respectively. Two genotypes at a QTL are the homozygote for the alleles inherited from Kefeng No.1 (solid) and the homozygote for the allele from Nannong 1138-2 (broken). Curves in grey are growth trajectories of 184 RILs. the kth trait were simulated by adding white noise with of noise variance σ were set as the estimates from the variances σ to the ODE curves for the jth QTL geno- 2 2 k real data, which are σ =2.42, σ =1.72 and 1 2 type with the probability w , i.e., the conditional prob- j|i 2 σ =0.14 for the leaf, stem, and root biomass, respec- ability of the ith RIL that carries the jth QTL genotype, tively. Meanwhile, by assuming a modest heritability given the two markers genotypes of this RIL. The values Biomass Biomass Wu et al. BMC Systems Biology 2011, 5:84 Page 6 of 11 http://www.biomedcentral.com/1752-0509/5/84 Kefeng No.1 Alleles Nannong 1138-2 Alleles Whole Whole Leaf Leaf Stem Stem Root Root 2 4 6 8 2 4 6 8 C D Whole Whole Leaf Leaf Stem Stem Root Root 2 4 6 8 2 4 6 8 Time Time Figure 3 Growth trajectories of whole-plant (red), leaf (green), stem (blue) and root biomass (black) for two different genotypes at a QTL detected on linkage group 12 (upper panel) and 24 (lower panel). Two genotypes at a QTL are the homozygote for the alleles inherited from Kefeng No.1 (A, C) and the homozygote for alleles inherited from Nannong 1138-2 (B, D). (0.05) for each trait at a middle stage of growth, we re- simulated data is a mimicry of the real soybean data, the results suggest that the experimental design used for scaled σ values which were used to simulate a new soybean mapping is scientifically sound and can provide data set. convincing QTL detection. This actually is not surpris- Systems mapping, implemented with the parameter ing because phenotyping has low measurement errors. cascading method, estimates QTL genotype-specific In analyzing a simulated data set for the traits curve parameters in the ODE (1) from the simulated assumed to have a modest heritability (0.05), the esti- data. The simulation was repeated 100 times to calculate mates of the ODE parameters are reasonably accurate the means, biases, standard deviations, and root mean and precise, indicating the power of systems mapping to square errors, with results tabulated in Table 2. It was detect small QTLs involved in trait formation. We per- found that the model can provide reasonably accurate formed an additional simulation to investigate the and precise estimates of QTL genotype-specific ODE power and false positive rates of the model by changing parameters with a modest sample size (n = 184). The levels of noises. In general, the power of the model is biases of the estimates are negligible, compared with the high, reaching 0.80 even when the heritability of growth scale of the standard deviations. Given that this Biomass Biomass Wu et al. BMC Systems Biology 2011, 5:84 Page 7 of 11 http://www.biomedcentral.com/1752-0509/5/84 Kefeng No.1 Alleles Nannong 1138-2 Alleles Leaf Stem Root Leaf Stem Root 2 4 6 8 2 4 6 8 C D Leaf Stem Root Leaf Stem Root 2 4 6 8 2 4 6 8 Time Time Figure 4 Time-dependent percentages of leaf (green), stem (blue) and root biomass (thin black) for two different genotypes at a QTL detected on linkage group 12 (upper panel) and 24 (lower panel). Two genotypes at a QTL are the homozygote for the alleles inherited from Kefeng No.1 (A, D) and the homozygote for alleles from Nannong 1138-2 (B, C). curves is low (0.05). Basically, a QTL can be fully statistical model called functional mapping has been detected when the heritability is 0.10 or larger. In any developed by incorporating the mathematical aspects of case, the false positive rates are not beyond 0.10, mostly trait development [19-21]. Despite significant improve- being less than 0.05. ment over conventional static mapping, functional map- ping has still a major limitation in characterizing Discussion developmental pathways that cause a final phenotype Mapping the genetic architecture of complex traits has and unraveling the underlying genetic mechanisms for been a subject of interest in both theoretical and empiri- trait formation and progression. cal aspects of modern biology [3-15]. Original In this article, we have for the first time put forward a approaches for genetic mapping are based on single- new approach-systems mapping by treating a phenotypic point variation in a phenotypic trait, neglecting the trait as a dynamic system and incorporating the design dynamic change of the trait during development. To principles of a biological system into a statistical frame- capture the dynamic pattern of genetic control, a new work for genetic mapping. Various components that Biomass Percentage Biomass Percentage Wu et al. BMC Systems Biology 2011, 5:84 Page 8 of 11 http://www.biomedcentral.com/1752-0509/5/84 Table 2 Means of maximum likelihood estimates of curve reveal the genetic control mechanisms for several parameters from the ODE system (1) and their biases, mechanistically meaningful relationships. They are (1) standard deviations (STD) and square root mean square size-shape relationship - is a big plant due to a big stem errors (RMSE) from 100 simulation replicates with sparse leaves or a small stem with dense leaves? (2) Genotype Estimate a b g a b a b g L L L S S R R R structural-functional relationship - in a specific environ- QQ TRUE 2.55 0.10 0.31 0.98 0.04 1.11 0.51 2.18 ment does a plant tend to allocate more carbon to its MEAN 2.55 0.10 0.31 0.98 0.04 1.11 0.51 2.18 leaves for CO uptake or roots for water and nutrient BIAS*10 2.10 -0.07 0.72 0.27 0.01 3.49 -0.88 -2.89 uptake? (3) cause-effect relationship - are more roots STD*10 3.75 0.08 0.60 0.94 0.01 1.29 0.49 1.11 duetomoreleavesordomoreleavesproducemore RMSE*10 3.75 0.08 0.61 0.94 0.01 1.34 0.50 1.15 roots? and (4) pleiotropic relationship - different traits with a similar function tend to integrate into modularity qq TRUE 1.89 0.14 0.44 1.04 0.07 1.11 0.56 1.85 [45]. How do the same QTLs pleiotropically control this MEAN 1.89 0.14 0.44 1.04 0.07 1.11 0.56 1.85 modularity? A better understanding of these relation- BIAS*10 -4.23 0.11 0.55 -0.61 -0.05 1.62 0.92 -2.02 ships helps to gain more insights into the mechanistic STD*10 3.60 0.14 0.92 0.99 0.04 1.38 0.45 1.66 response of plant size and shape to developmental and RMSE*10 3.62 0.14 0.92 0.99 0.04 1.39 0.46 1.67 environmental signals and, also, provide guidance to select an ideotype of crop cultivars with optimal shape and structure suited to a particular environment [46]. constitute the system through developmental regulation The model described in this article is a simple frame- are studied and connected by a system of biologically work for systems mapping. It can be used as a start meaningful differential equations (DE). Thus, the genetic mapping of a complex phenotype become an issue of point to expand the concept of systems mapping to testing and estimating genotype-specific curve para- tackle more complicated biological problems. A pheno- meters at specific QTLs that define the emerging prop- type can be dissected to any number of components at erties and dynamic behavior of the DE system. Systems any level of organization, molecule, cell, tissue, or whole mapping identifies QTLs that control developmental organism, depending on the interest of researchers and interactions of traits, the temporal pattern of QTLs data availability. With more knowledge about phenotype expression during development, as well as the genetic formation and development, more components can be determinants that control developmental switches (on/ involved in a system that is specified by high-dimension off). differential equations. Sophisticated mathematical tech- Systems mapping was applied to map QTLs for niques are needed to obtain stable solutions of these equations. In addition, by integrating it with genome- dynamic trajectories of biomass from different organs, wide association studies, systems mapping will not only the stem, leaves, and roots, that interact and coordinate provide a clear view of how different components inter- to determine whole-plant biomass growth, in an experi- act and coordinate to form a phenotype, but also will be mental cross of soybeans between Kefeng No. 1 and capable of illustrating a comprehensive picture of the Nannong 1138-2. In general, the stem receives a propor- genetic architecture of complex phenotypes. There is tionally larger amount of biomass with development, also a good reason to integrate systems mapping with accompanying a proportional decrease of biomass to the network biology to explore how “omics” information leaves and roots. Specific QTLs, biomass1 and biomass2, contribute to the regulatory mechanisms of phenotype that control this allometric change with development formation [47]. In any case, systems mapping will open have been detected from systems mapping. The alleles a new avenue for understanding the genetic architecture at the two QTLs inherited from Kefeng No. 1 tend to of complex phenotypes from a perspective of mechanis- amplify this contrast in development-dependent biomass partitioning, as compared to those from Nannong 1138- tic pathways inside their formation. 2. One of the two QTLs, biomass2, was found in a simi- lar genomic region identified by traditional functional Conclusions mapping [44]. This consistency does not only simply The past two decades have seen a phenomenal increase verify our systems mapping, but also gains new insight in the number of tools for the genetic mapping of com- into biological functions of the detected QTLs. For plex traits. Although genetic mapping continues to be example, the two QTLs detected, biomass1 and bio- an interesting area in genetic research owing to the suc- mass2, trigger genetic effects on the interactions and cess of molecular and sequencing technologies in gener- coordination of different organs which cause the ating a flood of data, a conceptual breakthrough in this dynamic variation of biomass growth. area remains elusive. In this article, we present a bot- Through various tests for ODE parameters individu- tom-top model for mapping and studying the genetic ally or in a combination, our systems mapping can architecture of complex traits. Different from existing Wu et al. BMC Systems Biology 2011, 5:84 Page 9 of 11 http://www.biomedcentral.com/1752-0509/5/84 mapping models, we use a systems approach to identify errors for such a large-scale field trial. Phenotyping pre- specific genes or quantitative trait loci that govern the cision was estimated to be above 95%. developmental interactions of various components com- Unlike a traditional mappingproject,ourgoal is to prising the phenotype. The map of developmental inter- map QTLs that control the dynamic process of how dif- actions among different components is constructed by a ferent organs, the stem, leaves, and roots, interact and system of differential equations. Thus, by estimating and coordinate to determine whole-plant biomass. The inter- testing mathematical parameters that specify the system, actions and coordination of different organs for a plant are understood using design principles described by the we are able to predict or alter the physiological status of ODE system (1). a phenotype based on the underlying genetic control mechanisms. We have tested and validated our model by analyzing a real data set for genetic mapping of bio- Likelihood mass growth in soybeans. The detection of QTLs by the Let denote the vector of y = y (t ) , ··· , y t ki i1 ki im ki i new model provides biologically meaningful interpreta- phenotypic values for trait k (k = 1 for leaf biomass (L), tions of QTL effects on trait formation and dynamics. 2 for stem biomass (S), and 3 for root biomass (R)) mea- The new model can be readily used to study the genetic sured on progeny i at time points .Note t , ··· , t i1 im basis of phenotypes in any other organism. i that the number of time points measured may be pro- geny-specific, expressed as m for progeny i.Assuming Methods that multiple QTLs (segregating with J genotypes), each Mapping Population bracketed by two flanking markers M, affects these Our model derivation is based on a mapping population three traits, we construct a mixture model-based likeli- comprising of n recombinant inbred lines (RILs), hood as initiated with two inbred lines. By continuous selfing or n J inbreeding, RILs after the F generation are considered L(z, M)= [ω f (z ;  , )] (2) homozygous because the fixation at any locus is given j|i j i j 7-1 i=1 j=1 by f = 1 - 0.5 ≈ 1. In practical terms, all plants from a single RIL are genetically identical, and can be used for where y =(y , y , y ) is a joint vector of phenotypic 1 2 3 replicated experiments under different environments. In values for the three traits, with z =(y , y , y )present- i 1i 2i 3i addition, each RIL represents a unique combination of ing the z-vector for progeny i; ω is the conditional j|i alleles from the parental genotypes where there are two probability of QTL genotype j (j = 1,..., J) given the mar- homozygous genotypes at each marker locus, each cor- ker genotype of progeny i, which can be expressed as a responding to a parental allele. The mapping population function of the recombination fractions between the is genotyped at molecular markers to construct a linkage QTL and markers [50], and f (z ; Θ , Ψ)isanMVN of j i j map covering the entire genome. The recombination leaf, stem and root biomass for progeny i which carries fraction between two markers is converted to the QTL genotype j, with mean vectors genetic distance in centiMorgan (cM) through a map function, such as the Haldane or Kosambi map function. μ =(μ , μ , μ ), j =1, ... , J j 1j 2j 3j The map constructed is used to locate QTLs that con- trol a quantitative trait of interest. specified by Θ , and covariance matrix specified by Ψ. We obtained a sample of 184 RILs derived from two If a system of differential equations (1) is used to jointly cultivars, Kefeng No. 1 and Nannong 1138-2, for map- model QTL genotype-specific means vectors for the ping agronomic traits. These RILs were genotyped for three traits, then we have Θ=(a , b , l , a , b , a , j Lj Lj jL Sj Sj Rj 950 molecular markers locating in 25 linkage groups b , l ) for genotype j. Rj Rj [48,49]. The plants were grown in a simple lattice design with two replicates in a plot at Jiangpu Soybean Experi- Estimation ment Station, Nanjing Agricultural University, China. Unlike a traditional mixture model for QTL mapping, Ten plants in the second row of a plot were randomly we will model the genotypic values of each QTL geno- type in likelihood (2) characterized by a group of non- selected for measuring leaf, stem and root biomass at linear ODEs. While an analytical solution is not each time in the whole growing season. After 20 days of available, we will implement numerical approaches to seedling emergence, dry weights separately for the solve these ODEs. Let μ (t) denote the genotypic value leaves, stem and roots were measured once every 5 to kj of the kth trait at time t for a QTL genotype j.Thus, 10 days until most plants stopped growth. A total of 6 the dynamic system of the traits and their interactions, to 8 measurements were taken for each of the RILs stu- died. Great efforts were made to control measurement regulated by QTL genotypes, can be modeled by a Wu et al. BMC Systems Biology 2011, 5:84 Page 10 of 11 http://www.biomedcentral.com/1752-0509/5/84 system of ODE (1), Additional material dμ (t) kj Additional file 1: Figure S1. The profiles of the log-likelihood ratios (LR) (3) = g (μ (t),  ), k =1, ··· ,3, k j between the full model (there is a QTL) and reduced model (there is no dt QTL) for soybean height growth trajectories throughout the soybean genome composed of 25 linkage groups. where μ (t)= (μ (t), ..., μ (t)) ,and Θ is a vector of j 1j 3j j ODE parameters associated with QTL genotype j.For J possible genotypes in the mapping population, we have T T Acknowledgements . The question now is how to Θ esti- =  , ... , 1 J This work is supported by the Changjiang Scholars Award, “One-thousand mate from noisy measurements. The functional mean Person” Award, NSF/IOS-0923975, a discovery grant of the Natural Sciences and Engineering Research Council of Canada (NSERC) (J. Cao), the National μ (t) may be represented as a linear combination of kj Key Basic Research Program of China (2009CB1184，2010CB1259， basis functions: 2011CB1093), The National High Technology R&D Program of China (2009AA1011), and the MOE 111 Project (B08025). Part of this work was carried when RW and JC were invited Research Fellows at the Statistical and μ (t)= c φ (t)= c φ (t) (4) kj kjr kjr Applied Mathematical Sciences Institute (SAMSI), sponsored by Duke kj kj University, University of North Carolina at Chapel Hill, and North Carolina r=1 State University. Thanks are due to Prof. Shouyi Chen and Prof. Deyue Yu for kind permission to use the jointly developed NJRIKY genetic linkage map. where j (t) = (j (t), ..., j (t)) is a vector of basis kj kj1 kjR functions with R orders and c =(c ,..., c ) is a vec- kj kj1 kjR Author details 3;J Center for Computational Biology, National Engineering Laboratory for Tree tor of basis coefficients. Define c = c as a kj k=1;j=1 Breeding, Key Laboratory of Genetics and Breeding in Forest Trees and Ornamental Plants, Beijing Forestry University, Beijing 100083, China. length (R ×3× J) vector of basis coefficients. The cubic Department of Statistics & Actuarial Science, Simon Fraser University, B-splines are often chosen as basis functions, since any Burnaby, B.C. Canada V5A 1S6. Department of Agronomy, Henan Institute of B-spline basis function is only positive over a short sub- Science and Technology, Xinxiang 453003, China. Center for Statistical Genetics, Pennsylvania State University, Hershey, PA 17033, USA. National interval and zero elsewhere. This is called the compact Center for Soybean Improvement, National Key Laboratory for Crop Genetics support property, and is essential for efficient computa- and Germplasm Enhancement, Soybean Research Institute, Nanjing tion. The flexibility of the B-spline basis functions Agricultural University, Nanjing 210095, China. Department of Horticultural Sciences, University of Florida, Gainesville, FL 32611, USA. depend on the number and location of knots we choose. It is an infinite-dimension optimization problem to Authors’ contributions choose the optimal number of knots and their locations. RW conceived of the idea of systems mapping, coordinated the whole study, and wrote the manuscript. JC derived the model, computed the real A popular approach to avoid this dilemma is choosing a data, run simulation studies, and participated in the writing of the Methods saturated number of knots and using a roughness pen- section. ZH conducted the soybean experiment and collected phenotypic alty to control the smoothness of the fitted curve and data. ZW packed the model into a computer package SysMap. JG directed the experimental design and data collection of soybeans. EV oversaw the avoid over-fitting [40]. project and highlighted the biological relevance a computational model We estimate the basis coefficient c and ODE para- must possess. All authors read and approved the final manuscript. meter Θ based on a two-nested level of optimization. In Received: 7 February 2011 Accepted: 27 May 2011 the inner level of optimization, c is estimated by opti- Published: 27 May 2011 mizing a criterion U(c|Θ), given any value of Θ.There- fore, the estimate maybeviewedasa function of Θ, cˆ References which is denoted as ˆ . Since no analytic formula 1. Lee SH, van der Werf JHJ, Hayes BJ, Goddard ME, Visscher PM: Predicting c () unobserved phenotypes for complex traits from whole-genome SNP for is available, is an implicit function. In the cˆ () cˆ data. PLoS Genet 2008, 4(10):e1000231.. outer level of optimization, Θ is estimated by optimizing 2. Lynch M, Walsh B: Genetics and Analysis of Quantitative Traits Sinauer Associates, Sunderland, MA; 1998. a criterion H , cˆ () . The parameter is removed in cˆ 3. Lander ES, Botstein D: Mapping Mendelian factors underlying quantitative traits using RFLP linkage maps. Genetics 1989, 121:185-199. the parameter space in the outer level by treating it as 4. Zeng ZB: Precision mapping of quantitative trait loci. Genetics 1994, an implicit function of Θ. Although cˆ () does not 136:1457-1468. have an analytic formula, the outer level of optimization 5. Xu S, Atchley WR: A random model approach to interval mapping of quantitative genes. Genetics 1995, 141:1189-1197. only requires to calculate the derivative dcˆ/d ,which 6. Wu RL, Ma CX, Casella G: Joint linkage and linkage disequilibrium can be obtained by using the implicit function theorem. mapping of quantitative trait loci in natural populations. Genetics 2002, 160:779-792. The above optimization procedure is called the para- 7. Yi N, Xu S: Bayesian lasso for quantitative trait loci mapping. Genetics meter cascading method. Note that when the two cri- 2008, 179:1045-1055. teria U(c|Θ)and H , cˆ () are the same, the 8. Zou F, Nie L, Wright FA, Sen PK: A robust QTL mapping procedure. J Stat Plann Infer 2009, 139:978-989. parameter cascading method is equivalent to the profil- 9. Ehrenreich IM, Torabi N, Jia Y, Kent J, Martis S, Shapiro JA, Gresham D, ing method. Caudy AA, Kruglyak L: Dissection of genetically complex traits with extremely large pools of yeast segregants. Nature 2010, 464:1039-1042. Wu et al. BMC Systems Biology 2011, 5:84 Page 11 of 11 http://www.biomedcentral.com/1752-0509/5/84 10. Steinmetz LM, Sinha H, Richards DR, Spiegelman JI, Oefner PJ, McCusker JH, 37. Huang JZ, Liu N, Pourahmadi M, Liu L: Covariance selection and Davis RW: Dissecting the complex architecture of a quantitative trait estimation via penalized normal likelihood. Biometrika 2006, 93:85-98. locus in yeast. Nature 2002, 416:326-330. 38. Li L, Brown MB, Lee KH, Gupta S: Estimation and inference for a spline- 11. Singer JB, Hill AE, Burrage LC, Olszens KR, Song J, Justice M, O’Brien WE, enhanced population pharmacokinetic model. Biometrics 2002, Conti DV, Witte JS, Lander ES, Nadeau JH: Genetic dissection of complex 58:601-611. traits with chromosome substitution strains of mice. Science 2004, 39. Ramsay JO: Principal differential analysis: Data reduction by differential 304:445-448. operators. J Roy Stat Soc Ser B 1996, 58:495-508. 12. Altshuler D, Daly MJ, Lander ES: Genetic mapping in human disease. 40. Ramsay JO, Silverman BW: Functional Data Analysis. 2 edition. Springer, New Science 2008, 322:881-888. York; 2005. 13. Shao H, Burrage LC, Sinasac DS, Hill AE, Ernest SR, O’Brien W, Courtland HW, 41. Ramsay JO, Hooker G, Campbell D, Cao JG: Parameter estimation for Jepsen KJ, Kirby A, Kulbokas EJ, Daly MJ, Broman KW, Lander ES, Nadeau JH: differential equations: a generalized smoothing approach (with Genetic architecture of complex traits: Large phenotypic effects and discussion). J Roy Stat Soc Ser B 2007, 69:741-796. pervasive epistasis. Proc Natl Acad Sci USA 2008, 105:19910-19914. 42. Liang H, Wu HL: Parameter estimation for differential equation models 14. Mackay TFC, Stone EA, Ayroles JF: The genetics of quantitative traits: using a frame-work of measurement error in regression model. JAm challenges and prospects. Nat Rev Genet 2009, 10:565-577. Stat Assoc 2008, 103:1570-1583. 15. Balasubramanian S, Schwartz C, Singh A, Warthmann N, Kim MC, Maloof JN, 43. Efron B, Tibshirani RJ: An Introduction to the Bootstrap Chapman & Hall; Loudet O, Trainer GT, Dabi T, Borevitz JO, Chory J, Weigel D: QTL mapping 1993. in new Arabidopsis thaliana advanced intercross-recombinant inbred 44. Li Q, Huang Z, Xu M, Wang C, Gai J, Huang Y, Pang X, Wu RL: Functional lines. PLoS ONE 2009, 4(2):e4318. mapping of genotype-environment interactions for soybean growth by 16. Frary A, Nesbitt TC, Grandillo S, Knaap E, Cong B, Liu J, Meller J, Elber R, a semiparametric approach. Plant Methods 2010, 6:13. Alpert KB, Tanksley SD: fw2. 2: a quantitative trait locus key to the 45. McCarthy MC, Enquist BJ: Consistency between an allometric approach evolution of tomato fruit size. Science 2000, 289:85-88. and optimal partitioning theory in global patterns of plant biomass 17. Li C, Zhou A, Sang T: Rice domestication by reducing shattering. Science allocation. Funct Ecol 2007, 21:713-720. 2006, 311:1936-1939. 46. Jiao Y, Wang Y, Xue D, Wang J, Yan M, Liu G, Dong G, Zeng D, Lu Z, Zhu X, 18. Huang X, Qian Q, Liu Z, Sun H, He S, Luo D, Xia G, Chu C, Li J, Fu X: Qian Q, Li J: Regulation of OsSPL14 by OsmiR156 defines ideal plant Natural variation at the DEP1 locus enhances grain yield in rice. Nat architecture in rice. Nat Genet 2010, 42:541-544. Genet 2009, 41:494-497. 47. Wu S, Yap JS, Li Y, Li Q, Fu GF, Li JH, Das K, Berg A, Zeng YR, Wu RL: 19. Ma CX, Casella G, Wu RL: Functional mapping of quantitative trait loci Network models for dissecting plant development by functional underlying the character process: a theoretical framework. Genetics 2002, mapping. Curr Bioinform 2009, 4:183-187. 161:1751-1762. 48. Li H, Huang Z, Gai J, Wu S, Zeng Y, Wu RL: A conceptual framework for 20. Wu RL, Lin M: Functional mapping - how to map and study the genetic mapping quantitative trait loci regulating ontogenetic allometry. PLoS architecture of dynamic complex traits. Nat Rev Genet 2006, 7:229-237. ONE 2007, 2(11):e1245. 21. Li Y, Wu RL: Functional mapping of growth and development. Biol Rev 49. Zhang W-K, Wang Y-J, Luo G-Z, Zhang J-S, He C-Y, Wu XL, Gai JY, Chen SY: 2010, 85:207-216. QTL mapping of ten agronomic traits on the soybean (Glycine max L. 22. Kitano H: Systems biology: a brief overview. Science 2002, 295:1662-1664. Merr.) genetic map and their association with EST markers. Theor Appl 23. Jansen RC: Studying complex biological systems using multifactorial Genet 2004, 108:1131-1139. perturbation. Nat Rev Genet 2003, 4:145-151. 50. Wu RL, Ma CX, Casella G: Statistical Genetics of Quantitative Traits: Linkage, 24. Csete ME, Doyle JC: Reverse engineering of biological complexity. Science Maps, and QTL Springer-Verlag, New York; 2007. 2002, 295:1664-1669. doi:10.1186/1752-0509-5-84 25. Cannell MGR, Dewar RC: Carbon allocation in trees: a review of concepts Cite this article as: Wu et al.: Systems mapping: how to improve the for modelling. Ad Ecol Res 1994, 25:59-104. genetic mapping of complex traits through design principles of 26. Luo Y, Field CB, Mooney HA: Predicting responses of photosynthesis and biological systems. BMC Systems Biology 2011 5:84. root fraction to elevated CO Interactions among carbon, nitrogen, and 2: growth. Plant Cell Environ 1994, 17:1195-1204. 27. Chen J, Reynolds J: A coordination model of carbon allocation in relation to water supply. Ann Bot 1997, 80:45-55. 28. Weiner J: Allocation, plasticity and allometry in plants. Perspect Plant Ecol Evol Syst 2004, 6:207-215. 29. Hermans C, Hammond JP, White PJ, Verbruggen N: How do plants respond to nutrient shortage by biomass allocation? Trends Plant Sci 2006, 11:610-617. 30. Marcelis LFM, Heuvelink E: Concepts of modelling carbon allocation among plant organs. In Functional-Structural Plant Modelling in Crop Production. Edited by: Vos J, Marcelis LFM, de Visser PHB, Struik PC and Evers JB. Springer. Printed in the Netherlands; 2007:103-111. 31. Genard M, Dauzat J, Franck N, Lescourret F, Moitrier N, Vaast P, Vercambre G: Carbon allocation in fruit trees: From theory to modelling. Trees - Structure and Function 2008, 22:269-282. 32. West GB, Brown JH, Enquist BJ: A general model for the origin of Submit your next manuscript to BioMed Central allometric scaling laws in biology. Science 1997, 276:122-126. and take full advantage of: 33. West GB, Brown JH, Enquist BJ: The fourth dimension of life: Fractal geometry and allometric scaling of organisms. Science 1999, • Convenient online submission 284:1677-1679. 34. West GB, Brown JH, Enquist BJ: A general model for ontogenetic growth. • Thorough peer review Nature 2001, 413:628-631. • No space constraints or color figure charges 35. Putter H, Heisterkamp SH, Lange JMA, De Wolf F: A Bayesian approach to • Immediate publication on acceptance parameter estimation in HIV dynamical models. Stat Med 2002, 21:2199-2214. • Inclusion in PubMed, CAS, Scopus and Google Scholar 36. Huang Y, Wu H: A Bayesian approach for estimating antiviral efficacy in • Research which is freely available for redistribution HIV dynamic models. J Appl Stat 2006, 33:155-174. Submit your manuscript at www.biomedcentral.com/submit

Journal

BMC Systems Biology – Springer Journals

Published: May 27, 2011

There are no references for this article.