Modeling the effects of N application on growth, yield and plant properties associated with the occurrence of chalky grains of rice

Abstract The objective of this study was to propose a model for explaining rice responses to a wide range of N application rates in various growth attributes associated with the occurrence of chalky grains. We improved the sub-model for N uptake process of a previous rice model which was originally developed for explaining genotypic and environmental variations in the whole growth processes, considering the difference in the rate of N loss from the plant-soil system between indigenously supplied soil mineral N and fertilizer N. A total of 80 growth datasets of cultivar ‘Koshihikari’ grown at Shiga prefecture, Japan, in 2010 was utilized for the calibration and validation of the model. The rice growth model well explained the above-ground biomass growth (RMSD = 78.7 g m−2) and rough dry grain yield (RMSD = 83.2 g m−2) for the validation data-set, simultaneously. The simulated carbohydrate content available per single spikelet was negatively correlated with the observed percentage of the milky-white grain which includes white-cored grain (r = −.77, p < .001) for all the data-sets of calibration and validation. On the other hand, the observed percentage of the sum of white-back and white-base grains was closely correlated with the simulated plant N content available per single spikelet (r = −.59, p < .001). It was suggested that the present rice growth model would rationally explain the effects of N application on the occurrence of the chalky grains through the dynamic change of the carbohydrate content and plant N content available per single spikelet.

The current global warming tendency (IPCC 2007) poses new challenges in rice production. It is well documented that the increase in air temperature above 26-27 °C during the 20-day period after the heading enhances the occurrence of chalky grains in rice (Morita et al., 2005). We need to develop an adaptive rice production technology to avoid the decline in rice quality due to the increased occurrence of chalky grains.
The percentage of chalky grain in rice was under the influence of not only high temperature but also plant properties such as spikelet number (Inoue, 2003;Tsukimori, 2003), the amount of carbohydrate available per single spikelet (Kobata et al., 2004), leaf color (leaf nitrogen (N) concentration) at around heading stage and N concentration of brown rice (Wakamatsu et al., 2008). A number of reports suggested that the short-term decline of carbohydrate supply to single grain due to the increased spikelet number and/or decreased carbohydrate content available for filling grains enhanced the occurrence of white-cored and milky-white grains (Kobata et al., 2004;Morita et al., 2005;Nakagawa et al., 2006;Nagato & Ebata, 1965;Tashiro & Wardlaw, 1991;Tsukaguchi & Iida, 2008). Furthermore, recent studies clearly showed that the percentages of white-back and white-base grains were more related to the plant N status than the spikelet number or the amount of carbohydrate available for filling grains; the decline in plant N concentration during the early stage of grain filling resulted in the enhanced occurrence of the white-back and white-base grains (Kusuda et al., 2004;Morita et al., 2005;Nakagawa et al., 2006;Tanaka et al., 2010).
All the above plant properties associated with the occurrence of chalky grains are determined by the dynamic change in the balance between carbon and N gains. It is well established that the spikelet number is closely correlated to the plant N content at around the late spikelet differentiation stage (Wada, 1969). Although the high leaf N content per unit leaf area is well accepted to increase the leaf photosynthetic rate, the high plant N concentration has also been associated with the declined non-structural carbohydrate concentration (Batten et al., 1993;Wada, 1981). Thus, the amount of carbohydrate available per single spikelet, defined as the sum of non-structural carbohydrate content accumulated up to the heading stage and biomass production during the period from heading to maturity divided by the number of spikelet (Kobata et al., 2004), is the final outcome of the dynamics in the balance OPEN ACCESS during the period from heading to maturity (Yoshida & Horie, 2009). The rate of NSC accumulation was described as a function of plant N concentration (Yoshida & Horie, 2009). The balance between sink (the product of spikelet number and single grain weight) and source sizes determined the growth rate of grains (Yoshida & Horie, 2009).
The model has been well validated for nine rice genotypes grown under seven climatically different sites in Asia (Yoshida & Horie, 2010). However, the explanatory power of the model for rice responses to various N application rates and timings has not been evaluated enough, because all of those rice growth data-sets utilized for calibrations and validations of the model were obtained only under the single protocol for fertilizer N application (Horie et al., 2003;Yoshida et al., 2006). To explain the rice responses to a wide range of N application rates in various growth attributes, we improved the sub-model for N uptake processes of the previous rice growth model as described below.

Improvement of the sub-model for N uptake processes
N uptake of rice under irrigated paddy field conditions was previously modeled considering the dynamics of soil mineral N, root system development toward horizontal direction and N-absorption activity of the root system (Yoshida & Horie, 2010). Although the previous model directly related the N uptake rate per unit root system to the amount of soil mineral N at that area, the model overestimated the growth response to the increased amount of basal and top-dressed fertilizer N for rice plants grown at Shiga prefecture in 2010 under various N application rates. In this study, the new model considered a difference in the rate of N loss from the plant-soil system between indigenously supplied soil mineral N and fertilizer N. Furthermore, the rate of plant N uptake was modeled based on the Michaelis-Menten kinetics between the rate of N uptake per unit root surface area and soil solution N concentration, the latter was estimated from the soil mineral N content per unit ground area considering the soil buffer power for ammonium N ion. Only these model modifications were proposed in this paper. For the complete structure of the sub-model for N uptake processes, please see Yoshida and Horie (2010).

Soil mineral N and soil solution N concentration.
The previous model assumed a one-layer root system of rice and described root system development toward horizontal direction as a function of lAI. To simulate the distribution of soil mineral N during root system expansion, the model divided the effective soil zone of a certain depth into 10 horizontal concentric partitions around an arbitrary rice hill, each partition having the same soil-water surface area. The new model divided the total soil mineral N content at ith soil partition (N soil(i) , between carbon and N gains. Therefore, the N fertilization method to optimize the spikelet number, nitrogen concentration and carbohydrate accumulation would be the central technology to avoid the decline in rice quality due to the increased occurrence of chalky grains. A crop growth model which explains the effects of N application rate and timing on growth, yield and plant properties associated with the occurrence of chalky grains in rice will be a useful tool to assess the adaptive rice production technology to a given environmental condition. The primary objective of this study is to propose a rice model for explaining rice responses to a wide range of fertilizer N application rate in various growth attributes associated with the occurrence of chalky grains. For this purpose, we improved the sub-model for N uptake processes of a previous rice model which was originally developed and proposed for explaining genotypic and environmental variations in the whole growth processes (Yoshida & Horie, 2010). A growth and yield data-set of rice cultivar 'Koshihikari' grown at Shiga prefecture Japan in 2010 was utilized for the calibration and validation of the improved sub-model for N uptake processes. Then, this paper shows the relationships between the observed occurrences of chalky grains and simulated variables associated with the occurrence of chalky grains for those rice plants grown at Shiga prefecture in 2010 under various N application rates. We discuss the possibility of predicting the occurrence of chalky grains based on a dynamic rice growth model.

Rice growth model
The process-based rice growth model proposed by Yoshida and Horie (2010) was used in this study. The details of the seven sub-models for explaining (1) phenological development, (2) plant N accumulation, (3) N allocation within plant organs, (4) leaf area index (lAI) development, (5) biomass growth, (6) spikelet number per unit area and (7) rough dry grain yield of irrigated rice, have already been presented in the previous papers (Yoshida & Horie, 2009Yoshida et al., 2006Yoshida et al., , 2007Yoshida et al., , 2008. Briefly, the model simulated the spikelet number per unit area as a function of plant N content two weeks before the heading and crop growth rate during the 2-week period preceding the heading; the former was associated with the number of differentiated spikelets and the latter was with the ratio of spikelet degeneration (Wada, 1969;Yoshida et al., 2006). The amount of carbohydrate that could be allocated to grains (source for yield formation) was predicted as the sum of non-structural carbohydrates (NSCs) accumulated in plants up to the heading stage and biomass production g m −2 ) into two N pools of N soilI (the pool of indigenously supplied soil N; g m −2 ) and N soilF (the pool of fertilized N; g m −2 ). Then, the daily change in soil mineral N at the ith soil partition from the hill (ΔN soil(i) , g m −2 day −1 ) was described as follows: where ΔN is(i) and ΔN f(i) represent the rates of indigenous N supply and N fertilization for the ith soil partition (g m −2 day −1 ); ΔN lossI(i) and ΔN upI(i) describe the rates of N loss and plant N uptake from the pool of N soilI at the ith soil partition (g m −2 day −1 ), and ΔN lossF(i) and ΔN upF(i) describe the rates of N loss and plant N uptake from the pool of N soilF at that soil partition (g m −2 day −1 ), respectively. The actual amount of N fertilization used in the field was entered in the model simulation and was equally distributed over 10 soil partitions as ΔN f(i) at each fertilization event. The rate of indigenous N supply was described as proportional to air temperature effective for N mineralization as described in Equation (3) of Yoshida and Horie (2010).
Because the rate of N loss due to ammonium volatilization, denitrification and leaching depends on the temperature, redox potential and other physical and chemical properties of the soil, the top-dressed N on the soil surface could be lost faster than the indigenously supplied N from the soil due to the environmental differences within soil layers (e.g. Reddy & Patrick, 1986;Schnier et al., 1988). Furthermore, because it takes a few days for fertilized N to be captured by the soil (De Datta et al., 1987), water management immediately after the N application could also enhance the N loss from the soil. Then, we assumed a difference in the rate of N loss from the plant-soil system between indigenously supplied soil N and fertilizer N. The ΔN lossI(i) and ΔN lossF(i) were most simply modeled as follows: where k SI and k SF are site-specific soil parameters representing the rate constant for N loss from indigenously supplied soil N and fertilizer N, respectively (day −1 ). The initial value of N soilI on the day of transplanting was estimated from the site-specific soil parameter for representing the rate of indigenous soil N supply at 25 °C (λ) as described in Yoshida and Horie (2010), while that of N soilF was the amount of basal fertilizer if any.
The dynamics of soil solution N concentration were estimated from that of soil N content per unit ground area, assuming 50 (cm 3 cm −3 ) for the soil NH 4 + buffer power ( Kirk and Kronzucker (2005)): where the ΔN soilconc(i) represents the daily change of soil solution N concentration at ith soil partition (μM), and constants 14 and 15 represent the atomic weight of N and the assumed soil depth (cm), respectively. The model ignored the gradients of soil solution N concentration within the rhizosphere.
N uptake rate per unit root system area. The previous model described the plant N uptake from ith soil partition from the hill (N up(i) , g m −2 d −1 ) as follows: where f(Δroot) represents the N-absorption activity of the root system which ranges from 0 to 1 and is described later and ω i the index of root system size which was defined as relative root system cover to the effective soil zone at ith soil partition. This Equation (7) in this paper can be obtained from the substitution of Equation (6) of Yoshida and Horie (2010). Briefly, our model adopted an index of root system size (ω) to describe root system development, which was modeled as a function of lAI utilizing a genotype-specific parameter α for representing the lAI at which root system fully covers effective soil zone (Equation (1) of Yoshida and Horie (2010): ω = lAI/α, 0 < ω ≤ 1). Setting the maximum ΔN up(i) to be 1.2 g m −2 day −1 (Peng & Cassman, 1998), this equation hypothesized that the N uptake rate per unit root system is proportionally increased with soil N content per unit ground area up to the maximum value. On the other hand, it has been well established that the relationship between solution N concentration at root surface and N uptake rate per unit root area is represented by the Michaelis-Menten kinetics (e.g. Youngdahl et al., 1982;Wang et al., 1993). The new model represented the ΔN up(i) as follows: where N maxup represents the maximum N uptake rate from the fully expanded root system (g m −2 d −1 ). The values of N maxup was estimated to be 1.71 by utilizing parameter values provided by Kirk and Kronzucker (2005) and where NSC H , ΔW HtoM and S represent the non-structural carbohydrate content accumulated up to the heading stage, biomass production during the period from heading to maturity and spikelet number per unit ground area, respectively, which were simulated by the rice growth model. Similarly, we defined the amount of N available per single spikelet (AN, mg) by the following equation: where NL H and LAI H represent the leaf N content (g m −2 ) and lAI at heading stage, respectively, and ΔN HtoM the N accumulation during the period from heading to maturity (g m −2 ), which were all simulated by the rice growth model. The constant 1.6 was utilized to estimate the total amount of available N from leaf N content, assuming that N was remobilized from leaves and stems simultaneously to grains in the ratio of 5:3 (see Equation (7) of Yoshida & Horie, 2009). The constant .45 represents the N content remaining in dead leaves (Yoshida & Horie, 2009).

Crop and weather data
Rice growth data from a field experiment conducted at Shiga Prefecture Agricultural Technology Promotion Center (136.13˚N) in 2010 were utilized in this study (Takehisa et al., unpublished). The field experiment was conducted on cultivar 'Koshihikari' for three transplanting dates (3 May, 17 May, and 27 May) with two replications. The effects of various N application rates and timings on the yield and occurrence of chalky grains were examined under two different soil types of clay loam (Cl) and sandy loam (Sl). As a total combination of three transplanting dates, two soil types and various N application rates and timings, 80 treatments were examined in the experiment (Table 1). A part of the results of this field experiment in 2010 has been reported in Hotta et al., (2011Hotta et al., ( , 2013 assuming 30 (cm cm −3 ) for the root length density of the fully expanded rice root system (Morita, 1988). The value of km (50 μM) was also quoted from Kirk and Kronzucker (2005). The previous model described the N-absorption activity of rice root systems, f(Δroot), as a function of air temperature and root growth rate (Δroot) per root system size (ω) (Equation (5) in Yoshid and Horie (2010)): where l(T) is the function to represent temperature effect on N absorption activity and b a genotype-independent parameter. The root growth rate was determined by canopy photosynthesis, maintenance respiration, the amount of sugar accumulated in above-ground biomass and the substrate allocation ratio to roots as a function of crop development index (DVI) (Yoshida & Horie, 2009). The DVI function for representing substrates allocation ratio to roots (f root (DVI)) was modified in this study from a linear function to an exponential decay function to fit it better to the data from Dingkuhn (1996): As a consequence, the value of parameter b involved in Equation (9) was re-calibrated as described later.
Finally, the rate of plant N uptake per unit ground area was given by integrating ΔN up(i) over 10 soil partitions.

Plant properties associated with the occurrences of chalky grains.
The amount of carbohydrate available per single spikelet for filling (AC, mg) was defined as the sum of non-structural carbohydrate content accumulated up to the heading stage and biomass production during the period from heading to maturity divided by the number of spikelet (Kobata et al., 2004;Tsukaguchi et al., 1996): a the basal and topdressing n were applied at around 3-5 days before and 20-31 days after the transplanting, respectively. b the fertilizer n of Panicle1 and Panicle2 were top-dressed at around 15-16 and 7-9 days before the heading, respectively.

Calibration and validation
The new plant N uptake process model contained five unknown parameters: a constant b for characterizing the response curve of root system activity, a genotype-specific parameter α for representing the lAI at which root system fully covers the effective soil zone (a genotype-specific parameter for 'Koshihikari'), a site-specific soil parameter λ for representing the rate of indigenous N supply at 25 °C (mg N m −2 day −1 ), a site-specific soil parameter k SI for the rate of N loss for the indigenous soil N (day −1 ) and a site specific soil parameter k SF for the rate of N loss from fertilized N (day −1 ). Out of the five parameters, the value of parameter λ for representing the rate of indigenous N supply was directly estimated to be 123.6 from the soil data on NH 4 -N content after 10-weeks incubation. The soil was sampled from a clay loam field where plants were transplanted on 3 May. Furthermore, the value of the genotype-specific parameter α for representing the lAI at which root system fully covers the effective soil zone was estimated to be .429 beforehand, based on the experimental data on the relationship between lAI development and soil solution ammonium N concentration for 'Koshihikari' grown at Tsukuba, Japan, in 2012, as the value of lAI when soil solution ammonium N at the center of four rice plants started to decline (Yoshida et al., unpublished). Values for all the rest three parameters of b, k SI and k SF involved in the N uptake process model were determined by the simplex method of estimation of parameters in a non-linear function (Haga & Hashimoto, 1980), so as to minimize the sum of square errors between measured and estimated dynamics of plant N accumulation. Of the total 80 plant data-sets, 54 data-sets of rice grown at clay loam fields were utilized for the model calibration. The other 26 data-sets of rice grown at sandy loam fields were utilized only for the model validation. The assumable variations in those site-specific parameters λ, k SI and k SF among soil types and transplanting dates were ignored in this study, because no significant difference was found in the soil NH 4 -N content after 10 weeks of incubation between those soil types and transplanting dates (data not shown).
was also applied at each fertilization event of basal and topdressing, while 1.25 times the amount of K 2 O to N was applied at the first (Panicle1) and the second (Panicle2) fertilization events during panicle formation period (see Table 1).
Plants were sampled at the stages of heading and maturity for determining the amount of above-ground biomass and plant N concentration. Furthermore, out of all the 80 treatments, additional measurements for the above-ground biomass and plant N concentration were conducted for 45, 58 and 57 treatments when 6th and 9th leaf blades on main stem were fully expanded and when panicle length reached 1 mm, respectively. Plant samples were oven dried at 80 °C for more than two days and their dry weights were determined. The plant N concentration was determined by the Kjeldahl method. Grain yield and yield components were determined at physiological maturity. The percentages of (1) milky-white grain including white-cored grain, (2) the sum of white-back grain plus white-belly grain, (3) white-base grain and (4) chalky grain as a total summing all (1) to (3) types of chalky grains were determined by a grain-quality testing machine (SATAKE, RGQI10B). In this study, the percentages of (1) milky-white grain including white-cored grain and (2) the sum of white-back plus white-belly grains were simply described as 'milky-white grain' percentage and 'white-back grain' percentage, respectively.
The model utilized the weather data of daily maximum and minimum air temperature and daily total solar radiation to simulate the whole rice growth processes. All of those weather data required for the model simulation were obtained from the local meteorological observatory located at Hikone in Shiga prefecture through the MeteoCrop DB developed by NIAES (2009). Importantly, the observed daily mean air temperature in August and September was markedly higher than in normal year at the site Hikone in 2010. The monthly averages of daily mean air temperature on August and September were 29.3 °C and 25.0 °C in 2010, respectively, while those values of normal year were 27.1 °C and 23.2 °C, respectively. growth attributes: the dynamics of plant N accumulation, above-ground biomass growth, plant N concentration, and spikelet number per unit ground area and rough dry grain yield. The goodness of fit between measurements and simulations for each of the five objective variables was examined by means of RMSD, root mean squared deviation, and relative RMSD (rRMSD), the ratio of RMSD to mean of all the observations. Furthermore, we calculated the slope and coefficient of determination (R 2 ) of the regression line forced to go through the origin and correlation coefficient (r) between the measurements and simulations for each objective variable. Figure 1 shows the relationship between the measured and simulated dynamics of plant N accumulation in total above-ground biomass for (A) calibration data-set and (B) validation data-set, respectively. The new sub-model for N uptake process well explained the dynamics of plant N accumulation for 'Koshihikari' grown at Shiga prefecture in 2010 under various N application rates. The determination coefficient (R 2 ) and slope of the regression lines forced to go through the origin were .98 and 1.00 for the calibration data-set and .95 and 1.02 for the validation data-set, respectively. The estimated values of soil parameters k SI which represent the rate of N loss from the soil-plant system for indigenous soil N and k SF for fertilizer N were .026 and .074, respectively. Table 2 summarizes the validity of the model for five objective variables. The rice growth model integrated with the new sub-model for N uptake process well explained not only the dynamics of plant N accumulation but also the dynamics of the above-ground biomass growth and plant N concentration, simultaneously. The RMSD and rRMSD values between the measurements and simulations were 80.5 g m −2 and .13 for the above-ground biomass growth and .29 (g g −1 , %) and .17 for the plant N concentration for the calibration data-set, respectively (Table 2). However, the model underestimated the spikelet number per unit area and hence the rough dry grain yield; the slope of regression line forced to go through the origin was .84 for the spikelet number and .76 for the rough dry grain yield for the calibration data-set, respectively. These results shown in Table 2 were obtained under an assumption that values of genotype-specific parameters of 'Koshihikari' were all the same as those of cultivar 'Nipponbare' , except for those associated with the phenological development and the root system expansion (α).

Model performance
A comparison of Tables 2 and 3 shows the validity of the model for the five objective variables obtained after the The rice growth model contained 10 genotype-specific parameters to describe the genotypic differences in growth and yield formation processes, besides one parameter α for N uptake processes (Yoshida & Horie, 2009). In this study, values of genotype-specific parameters for 'Koshihikari' involved in the sub-model for phenological development were quoted from Horie et al., (1995), and those in the other sub-models, such as stomatal conductance, light extinction coefficient in the canopy and single grain weight, were all approximated to be the same as those values of cultivar 'Nipponbare' . However, because the model underestimated the spikelet number per unit ground area (see Table 2 at result section), we calibrated the value of the genotype-specific parameter A of 'Koshihikari' which represents N use efficiency for spikelet formation (Yoshida et al., 2006). This calibration was made by simplex method so as to minimize the sum of square errors between measured and simulated spikelet number per unit ground area, utilizing the same 54 plant data-sets of model calibration.
The evaluation of the rice growth model was conducted by comparing the measured and simulated values of five   three site-specific soil parameters of λ for representing the rate of indigenous N supply at 25 °C, k SI for the rate of N loss for the indigenous soil N and k SF for the rate of N loss from fertilized N (Table 3, Figures 2 and 3). Figure 4 shows the relationships between the simulated carbohydrate content available per single spikelet (AC, simulated by Equation (11)) and the observed milkywhite grain percentage for (A) calibration data-set and (B) validation data-set, respectively. Those simulations were conducted after model calibrations (see Table 3). The simulated carbohydrate content available per single spikelet was negatively correlated with the observed milky-white grain percentage: r = −.77, p < .001, for calibration data-set and r = −.81, p < .001, for validation data-set, respectively. Table 4 summarizes the correlation coefficients between the observed occurrences of chalky grains and simulated growth attributes, which contains both data-sets utilized for the model calibration and validation. The occurrence of milky-white grain was most closely correlated with AC (r = −.77, p < .001), followed by the plant N content available per single spikelet (AN, simulated by Equation (12); r = −.75, p < .001) and spikelet number per unit area (r = −.63, p < .001). On the other hand, the percentages of white-back and white-base grains were more closely correlated with AN (r = −.62, p < .001, for white-back grain and r = −.57, p < .001, for white-base grain, respectively) than AC (r = −.56, p < .001, for white-back grain and r = −.32, p < .01, for white-base grain, respectively). The percentage of chalky grains as a total was most closely correlated with the plant N content available per single spikelet with the correlation coefficient of −.72, p < .001 (Table 4). Figure  5 shows the relationship between the simulated plant N content available per single spikelet (AN) and the observed percentage of chalky grains as a total for (A) calibration data-set and (B) validation data-set, respectively.

Discussion
To simulate the root system development as a function of lAI, the model utilized a genotype-specific parameter α for representing the lAI at which root system fully covers the effective soil zone. Although we estimated the α value of 'Koshihikari' to be .429 based on the experimental data on the relationship between lAI development and soil solution ammonium N concentration (Yoshida et al., unpublished), it was far smaller than those values estimated by the previous N uptake sub-model for nine rice genotypes which ranged from 1.82 to 2.73 (Yoshida & Horie, 2010). In this study, the model predicted that it took 5-6 weeks for calibration of one genotype-specific parameter A for representing the N use efficiency for the spikelet production. By the calibration of the parameter A, the validity of the model for explaining the spikelet number per unit area and the rough dry grain yield was improved with RMSD and rRMSD values of 1971 (m −2 ) and .07 for the spikelet number per area and 70.1 (g m −2 ) and .11 for the rough dry grain yield for the calibration data-set, respectively (Table 3). Figures 2 and 3 show the relationships between the measured and simulated (A) dynamics of the aboveground biomass growth, (B) plant N concentration, and (C) spikelet number per unit area and (D) rough dry grain yield for the calibration data-set and validation data-set, respectively. Although the R 2 value of the regression line forced to go through the origin was poor for the spikelet number per unit area for the validation data-set (R 2 = .24), the rRMSD value .07 (smaller than 10%) was considered to be satisfactorily low (Jamieson et al., 1991). The rice growth model integrated with the new sub-model for N uptake process well explained the rice responses in those growth attributes to various N application rates and timings for cultivar 'Koshihikari' grown at Shiga prefecture, utilizing by the previous report that the soil solution ammonium N concentration decreased and became almost zero at around 40 days after the transplanting or the plant stage of mid-tillering at the center point of surrounding four plant hills, presumably reflecting the completed expansion of the rice root system to the farthest point from the arbitrary rice hill at that stage (Toriyama & Ishida, 1987). The difference in the parameter α among genotypes, lower value of Koshihikari in this study than those of other genotypes in the previous report (Yoshida & Horie, 2010), might be derived from the structural difference in the sub-model for N uptake processes. To compare the characteristics of 'Koshihikari' by utilizing the genotype-specific parameter α with other genotypes, we need further validations and calibrations of the model integrated with the new N uptake sub-model against various rice genotypes grown under wide range of soil N conditions.
From the model calibration, we estimated a value .026 for the site-specific soil parameter k SI for representing the rate of N loss from the indigenous soil N pool and .074 for the other soil-specific parameter k SF for representing that from the fertilizer N pool for the experimental field in Shiga prefecture, Japan. As expected, the value of k SF was larger than that of k SI . As a consequence, the simulated crop recovery of fertilizer N, which is comparable to the experimental data with 15 N labeled fertilizer was ranged from 26-48% for those rice plants grown at Shiga prefecture under various N application rates and timings. Those values were within the range or smaller than the previously reported fertilizer N recovery by various rice genotypes grown under wide range of environmental conditions (30-80% : Cassman et. al., 1993 and2002;Peng et. al., 1998, Dobermann et al., 2000. The range of the crop recovery of fertilizer N simulated for our experiments seemed to be realistic because the slopes of regression lines between the amount of total N fertilizer and plant N content at maturity were comparable between the observation (.42) and simulation (.45), combining all the data-sets of calibration and validation ( Figure 6). Figure 6 showed the observed and simulated plant N content at maturity under given N fertilization rates. Furthermore, the new N uptake process model largely improved the explanatory power for the rice responses to N application comparing with the original model with unimproved N uptake process model (Yoshida & Horie, 2010), which underestimated the rice N accumulation under low N application rate and overestimated it under high N application (the slope of the regression line between the amount of total N fertilizer and plant N content at maturity was 1.22 for the original model). Although this study ignored the possible differences in site-specific soil parameters between soil types of Cl and Sl, the parameter values of k SI and k SF can be different between the rice root system to cover the effective soil zone (data omitted), utilizing the value .429 for the genotype-specific parameter α. This simulation result is well supported  properties on the value of site-specific soil parameters k SI , k SF and λ are the subjects for further research. soil types, respectively, due to the large difference in soil properties even at the same location. The effects of soil Figure 7. the relationship between the simulated carbohydrate content available per single spikelet and the amount of fertilizer n applied at (a) basal: 3-5 days before transplanting, (B) topdressing: 20-31 days after the transplanting, (c) panicle1: 15-16 days before the heading, and (d) panicle2: 7-9 days before the heading for all the data-sets of calibration and validation. the figure legends show the transplanting date. Table 3. the goodness of fit between measurements and simulations for each variable after the calibration of a genotype-specific parameter. a rMSd and rrMSd represent the root mean squared deviation between measurements and simulations and relative rMSd, respectively.  (Table 4). This simulation results might be supported by recent reports that the percentage of whiteback and white-base grains was more related to the plant N status than the spikelet number or the amount of carbohydrate available for filling grains (Kusuda et al., 2004;Morita et al., 2005;Nakagawa et al., 2006). utilizing the above two plant factors of carbohydrates and plant N content available per single spikelet, the rice growth model may be able to rationally explain the occurrences of chalky gains of rice grown under various N application rates. However, it should be noted that the most critical factor that triggers the occurrence of chalky grains would be the high air temperature during the grain filling period (Morita et al., 2005). As mentioned before, the observed daily mean air temperature in August and September was markedly higher than in normal year at Shiga prefecture in 2010. Although the field experimental data utilized in this study for the model calibration and validation contained a variation in transplanting date, the observed air temperature during the 20-day period after the flowering was not differed among those rice growing seasons and notably high, ranging from 28.53 (for rice plants transplanted on 3 May), 28.96 (for those transplanted on 17 May) to 29.01 °C (for those transplanted on 27 May). We found that plant factors of carbohydrates and plant N content available per single spikelet would be good intermediate variables for the rice growth model to explain the effects of N application rates and timings on the occurrence of chalky grains. However, further research is needed to develop a more comprehensive model for explaining the occurrence of chalky grains considering the effects of air temperature during the grain filling stage.

Acknowledgment
I would like to thank Dr. Xinyou Yin, senior scientist of Centre for Crop Systems Analysis, Wageningen university who gave me constructive comments and Mr. Satoru Hotta who conducted the field experiment.

Disclosure statement
No potential conflict of interest was reported by the authors.

Funding
This work was financially supported by the Ministry of Agriculture, Forestry and Fisheries, Japan through a research project entitled "Development of technologies for mitigation and adaptation to climate change in Agriculture, Forestry and Fisheries". A part of this work was also supported by MEXT KAKENHI [grant number 22580021]. Figure 7 shows the relationship between the simulated carbohydrate content available per single spikelet and the amount of fertilizer N applied at (A) basal: 3-5 days before transplanting, (B) topdressing: 20-31 days after the transplanting, (C) panicle1: 15-16 days before the heading, and (D) panicle2: 7-9 days before the heading for all the data-sets of calibration and validation. The simulated carbohydrate content available per single spikelet was negatively correlated with the amount of fertilizer N at panicle1 (r = −.88, p < .001), while it was positively correlated with that at panicle2 (r = .82, p < .001). The decreased carbohydrate content per single spikelet with increase in fertilizer N at panicle1 was explained by the positive response of the spikelet number per unit area to the increase in fertilizer N at panicle1; the model predicted the number of differentiated spikelet as a function of plant N content 2 weeks before the heading, in between the dates of panicle1 and panicle2 (Yoshida et al., 2006). It should be noted that the positive correlation between the carbohydrate content available per single spikelet and the amount of N fertilizer at panicle2 was mainly derived from the interdependency between the amounts of fertilizer N at panicle1 and panicle2 (Table  1). A number of reports suggested that the short-term decline of carbohydrate supply to single grain due to the increased spikelet number and/or decreased carbohydrate content available for filling grains enhanced the occurrence of white-cored and milky-white grains (Kobata et al., 2004;Morita et al., 2005;Nagato & Ebata, 1965;Nakagawa et al., 2006;Tashiro & Wardlaw, 1991;Tsukaguchi & Iida, 2008). The model also predicted a close correlation between the carbohydrate content available per single spikelet and the occurrence of milky-white grain ( Figure 4, Table 4), which is well supported by those previous reports. Furthermore, the observed occurrence of milky-white grain was more closely correlated with the simulated carbohydrate content available per single spikelet (r = −.77, p < .001) than simulated (r = .63, p < .001) spikelet number per unit area. It was suggested that the present rice growth model integrated with the new N uptake sub-model would rationally explain the effects of N application rates and timings on the occurrence of the milky-white grain through the dynamic change of the carbohydrate content available per single spikelet.
This simulation study also suggested that the percentage of the sum of white-back and white-base grains was closely correlated with the plant N content available per single spikelet (r = −.59, p < .001), followed by the plant N concentration at maturity (r = −.39, p < .001), the amount of carbohydrate available per single spikelet (r = −.34, p < .01) and the spikelet number per unit area (r = .22,