Modeling of radial growth curves and radial variation of basic density in Chamaecyparis obtusa planted in two progeny test sites

The objectives of the present study are to clarify the effect of macro-and micro-environment on the radial growth patterns and radial variation patterns of basic density in hinoki cypress ( Chamaecyparis obtusa (Sieb. et Zucc.) Endl.). We evaluated the radial variation patterns of cumulative annual ring width (as radial growth pattern) and basic density by modeling methods using hinoki cypress 36 families planted at two progeny test sites. In addition, narrow-sense heritability and correlation between sites for annual ring width and basic density were investigated. As the results of modeling for radial growth patterns, radial growth patterns slightly differed between sites. In addition, the stem diameter reaching the plateau might be varied among blocks in a site. On the other hand, radial variation of basic density was affected by genetic factors rather than blocks in the site. However, the radial growth rate may somewhat affect the radial variation of basic density. The heritability and correlation coefficients between sites in basic density were higher than those of annual ring width. Therefore, although radial growth in hinoki cypress varies by the effects of micro-and macro-environmental factors and has some influence on the radial variation of basic density, basic density is more strongly affected by genetic factors than by these influences, allowing for effective improvement for wood density by tree breeding program.


Introduction
Wood density is one of the essential wood properties because wood density is positively correlated with the mechanical properties of wood [1][2][3][4][5].Wood properties, including wood density, vary from pith to bark, known as radial variation [6].If wood properties dramatically change from pith to bark, the lumber properties obtained from a log may vary.Therefore, producing a log with more uniform wood properties from pith to bark is desirable in the wood industry.However, the radial variation pattern may differ among genetic backgrounds and environmental factors within a species [2,5,7,8].If larger genetic effects on radial variation patterns of wood properties are found, tree breeding will be able to reduce the variation from pith to bark.
The wood near the pith where wood properties are more variable is called juvenile wood, while the outer side, which exhibits relatively stable wood properties, is defined as mature wood.The boundary cambial age between juvenile and mature wood in softwood species had been decided by radial variation of tracheid length [5,9,10] and/or other wood properties such as wood density [7,11].Whether or not the boundary cambial age between juvenile and mature wood varied by the effects of genetic factors and environmental factors in the softwood species should be investigated.
The growth characteristics and wood properties of trees were affected by genetic factors, environmental factors, and their interactions [12].The larger influence of genetic factors on the target traits, the greater potential for improvement by the tree breeding programs [13].For environmental factors, there are the macro-environment, such as climatic conditions, and the micro-environment, such as slope direction and inclination, soil type, and soil water content.In addition, interactions between genetic and environmental factors, known as Genotype (G) × Environment (E) interaction, could be detected by comparing the performances of genotypes among different environments [14].If there is a substantial effect of G × E interaction, growth characteristics and wood properties of genotypes might vary among different environmental conditions.Therefore, understanding the impact of the genetic, environmental, and interactional effects on target traits is vital to promoting efficient tree breeding programs.
Growth curves are expressed in various sigmoidlike models such as logistic, Chapman-Richards, and Gompertz functions, and they help predict the stem size at certain ages and characterize growth patterns [15][16][17].Zeltiņš et al. [17] pointed out that growth curve patterns varied among clones in Picea abies.Based on the asymptote and rate parameters of growth curves obtained from the logistic function applied to growth curves, Nagamitsu et al. [15] found that Larix kaempferi trees originated from different provenances exhibited the stable radial growth patterns across different provenance test sites in Japan.To increase wood volume production effectively, it is important to clarify whether the differences in genotypic, environmental and their interaction factors affect the growth curve patterns of stem in the target species.
Chamaecyparis is a genus of the Cupressaceae family distributed in the coastal areas of the Pacific and Atlantic Oceans in the northern hemisphere and is native to eastern Asia and the western and eastern margins of North America [18].The Chamaecyparis obtusa (Sieb.Et Zucc.)Endl.(hinoki in Japanese and hinoki cypress in English) is one of the most important forestry species in Japan, accounting for 25% of the total plantation area.The wood of this species has been used for construction materials [19], leading to the fact that wood properties and growth characteristics are important in this species.In the tree breeding program conducted in Japan, the evaluation of growth characteristics and wood properties of the second-generation plus tree clones originated from the first-generation plus tree clones have been proceeded [2,5,20].In the present study, we collected 985 core samples from 36 hinoki cypress families planted at two progeny test sites.To evaluate the radial variation patterns of cumulative annual ring width (as radial growth pattern) and basic density, we then analyzed their radial variation patterns by a nonlinear mixed-effects modeling approach.The objectives of the present study are to clarify the effect of macro-and micro-environment on the radial growth patterns and radial variation patterns of basic density.In addition, we investigated the radial variation of genetic control for annual ring width and basic density.

Progeny test sites
The hinoki cypress (Chamaecyparis obtusa (Sieb.Et Zucc.)Endl.) trees assessed in the present study were sampled from two progeny test sites located in Hitachi city, Ibaraki prefecture (Hitachi site) and Kiso district, Nagano prefecture (Kiso site) in Japan (Fig. 1).Hitachi site had warmer and drier conditions, whereas Kiso site had a colder, more rainy in summer, and snowy in winter (Table 1).In addition, the Kiso site is located at a high elevation with a relatively steep slope (Table 1).Trees in the Kiso site suffered from substantial animal damage by deer and meteorological damage, bears, and caused basal stem bending due to snow pressure.
Two progeny test sites have consisted of 36 common full-sib and open-pollinated families originating from Fig. 1 Location of the two-progeny test sites and locations where the plus trees used in the present study were selected.The black triangle locates the progeny test sites.The grey circle indicates the location where plus trees used as mating parents in the present study were selected 33 first-generation plus-tree clones.First-generation plus-tree clones used for those families in the present study were selected from various places in the Kanto and Chubu regions in Japan (Fig. 1).Hitachi and Kiso sites were established in April and May in 1998, respectively, with 3000 stems/ha as initial planting density, using a randomized block design with six replicates (Table 1).In the present study, 25 full-sib families and 11 open-pollinated families with sufficient individuals were subjected.In the Hitachi site, all trees were subjected, while around ten trees without severe damages due to animals and snow were selected per a family in the Kiso site.Finally, the total number of target individuals at the Hitachi and Kiso sites were 595 and 390 trees, respectively.

Annual ring width and basic density
Core samples of 12 mm in diameter from bark to pith were collected at 1.3 m above the ground from all the individuals using an increment borer (Haglöf Sweden, Långsele, Sweden).The cross sections of the obtained cores were planed by a core microtome (WSL, Birmensdorf, Switzerland).After that, cross sections were scanned with 1200 dpi resolution using a scanner (GT-9300UF, EPSON, Nagano, Japan).Annual ring width was measured using the cross-sectional images by the image analysis software (ImageJ, National Institute of Health, Bethesda, MD).The boundaries of annual rings were determined visually by the color shade of the latewood.The cumulative annual ring width (CRW) of each annual ring was calculated, which is used for modeling the growth curve described below.All core samples were cut at every three annual rings from the pith to bark to determine the basic density.Basic density was determined by dividing oven-dry weight at 105 ℃ by green volume determined by water displacement.Area-weighted basic density in an individual was calculated using the annual ring width dataset.Differences between sites were assayed for stem diameter, annual ring width mean from pith to bark, and area-weighted basic density by t-test.

Nonlinear modeling
Statistical analysis was performed using R version 4.1.2[22].To model radial growth pattern and radial variation of basic density in relation to annual ring number, nonlinear mixed-effect models were developed using the nlme function in the nlme package [23] (Table 2).To evaluate the radial growth patterns, cumulative annual ring width in relation to annual ring number from pith was first fitted using nonlinear models based on the Chapman-Richards function (Model I in Table 2).In this model, β 0 , β 1 , and β 2 are the fixed-effect parameters for the asymptotic diameter parameter, the rate parameter, and the shape parameter, respectively (see Fig. 2 for the response due to changes in each parameter) [17].Radial variation of basic density in relation to annual ring number was fitted using the nonlinear model based on the quadratic function (Model II in Table 2).In this model, β 0 and β 1 are the fixed-effect parameters for shape and variability of radial variation, respectively.In contrast, the β 2 is an intercept parameter (Fig. 2).Next, the nonlinear models were

Model Formula
then updated to include the random-effect parameters accounting for the variation of family or block into each model term to evaluate the effects of genetic or microenvironmental variation on the radial growth patterns and basic density (Table 2).All nonlinear mixed-effect models were separately fitted for each site to improve model convergence for robust parameter estimates.The model showing the lowest Akaike information criterion (AIC) [24] was selected as the best model.Furthermore, the annual ring number from pith reaching 1% of the variability in each site was calculated as the boundary between juvenile and mature wood using the fixed-effect parameters of the best model.

Genetic statistical analysis
Genetic variance for each property by sites was estimated with a linear mixed-effects model with restricted maximum likelihood estimation (REML, method = "em") using remlf90 function of breedR package [25] in R following formula: (1) where y ijkl is the measured value of the lth individual of the jth and kth parents in the ith block, µ is the general mean, B i is the fixed effect of the ith block, A jkl is the random effect of the additive genetic effect of the lth individual of the jth and kth parents, F jk is the random effect of family of the jth and kth parents, i.e., specific combining ability (SCA), and e ijkl is the residual.The random factors and breeding values were obtained using an individualtree model of the best linear unbiased prediction (BLUP).
The narrow-sense heritability (h 2 ) of annual ring width and basic density at each radial position was estimated using the following formula: where σ a 2 , σ f 2 , and σ e 2 are the variance of additive genetic effect, the variance of family and residual, respectively.The genetic variance, breeding value, and narrow-sense heritability were estimated for each site.
To evaluate the effect of interaction between genetics and environment on all properties and radial variation patterns of interaction, the correlation coefficient between sites was calculated using mating-parental (2) breeding values or random effects of family in the best models.Low correlation coefficients between sites imply a greater influence of the interaction of genetic and macro-environment [14].

Mean values in each site
The mean stem diameter in both sites was almost similar (15.7 cm in Hitachi and 15.6 cm in Kiso) (Table 3).
On the other hand, the annual ring width at the Kiso site (3.43 mm) was larger than that of the Hitachi site (3.19 mm).Area-weighted basic density at Kiso (0.389 g/ cm 3 ) was also higher value than that of Hitachi (0.365 g/ cm 3 ).The mean values of annual ring width and basic density obtained in the present study were consistent with previous reports [2,[26][27][28].The mean of annual ring width and area-weighted basic density significantly differed between sites.

Modeling of radial growth and radial variation of basic density
The AIC values of all developed models are shown in Table 4.The model I-i-b was selected as the best model for radial growth in both sites.This model had a random effect of block in the asymptotic diameter parameter term.As to the radial variation of basic density, model II-v-f was selected as the best model in both sites based on the AIC value.This best model included the random effect of the family in the variability and intercept parameter terms.Fixed-effect parameters of the selected models are listed in Table 5.All estimated parameters were significant.In the models for radial growth, asymptotic (β 0 ), rate (β 1 ), and shape (β 2 ) parameters were 78.6, 0.096, and 1.293 for Hitachi, and 85.7, 0.093, and 1.541 for Kiso, respectively.On the other hand, the shape (β 0 ) and variability (β 1 ) parameters of models for basic density for Hitachi were 0.000381 and −0.0133, similarly for Kiso, 0.000295 and −0.0108, respectively.The intercept parameters (β 2 ) for Hitachi and Kiso were 0.464 and 0.475, respectively.Based on the estimated parameter values in the best model, radial growth curves and radial variation of basic density are illustrated in Fig. 3.The radial variation pattern of basic density observed in the present study was almost consistent with previous reports in hinoki cypress [5,29,30].Annual ring numbers from the pith exhibiting an annual changing ratio of basic density reaching less than 1% were 13.9 at Hitachi and 12.5 at Kiso.

Heritability and correlation coefficients between sites
The narrow-sense heritability of annual ring width and basic density were estimated at each radial position (Table 6).Heritability of annual ring width at Hitachi   5) and random effects of family or block in sites for the best-fitted models selected by AIC (Table 4).The vertical solid lines indicate the annual ring number from the pith in which the rate of change in basic density reached less than 1% (13.9 years at Hitachi and 12.5 years at Kiso) increased from the pith side (h 2 = 0.003) to the bark side (h 2 = 0.386), whereas that of Kiso showed relatively low values (from h 2 = 0.005 to h 2 = 0.119) in almost all radial positions.On the other hand, the heritability of basic density in both sites was higher values at almost all radial positions (h 2 : from 0.094 to 0.502 at Hitachi; h 2 : from 0.188 to 0.465 at Kiso) compared to those of annual ring width.Similar trends were also found in correlations between sites: correlations between sites of basic density were higher at almost all radial positions than those of annual ring width in the range from 0.477 to 0.632.The correlation coefficients for b 1 in selected models for basic density between sites were not significant (r = 0.256), while that for b 2 is significant (r = 0.653) (Fig. 4).

Discussion
In the present study, our nonlinear modeling approach successfully modeled the radial patterns in the cumulative ring width (radial growth curve) and basic density, and quantified the genetic and environmental effects on those patterns.No significant difference between sites was found in stem diameter, while the mean annual ring width at the Kiso site was significantly higher than that of the Hitachi site (Table 3), suggesting that radial growth patterns might differ between sites.In the present study, cumulative annual ring widths were modeled in relation to annual ring numbers from pith to assess the difference between sites of radial growth patterns.As a result, the radial growth pattern slightly differed between sites.For example, in the comparison between sites of radial growth patterns based on the fixed-effect parameters, shape parameter (β 2 ) for Hitachi was lower than that for Kiso (Table 5), indicating that the radial growth rate at younger age at Hitachi was more rapid than that of Kiso (Figs. 2 and 3).However, the growth at Kiso was sustained for a longer time than at Hitachi since the rate parameter (β 1 ) for Kiso was lower than that for Hitachi.Therefore, the radial growth patterns of the hinoki cypress might be slightly affected by macro-environmental factors such as climatic conditions.In addition, the model included a random effect of progeny test's block in the asymptote term (Model I-i-b) was selected at both sites (Table 4), suggesting that the stem diameter reaching the plateau might be varied among blocks of each site.Nishizono et al. [31] reported that the strong influence of temperature on radial growth in C. japonica planted at sixteen sites throughout Japan.On the other hand, Zacharias et al. [32] mentioned that micro-environment such as nutrient availability, water holding capacity of the soil, and light conditions were likely to be more relevant for tree growth than genetic similarity in Picea glauca.In the hinoki cypress, the radial growth rate might be affected by macro-and micro-environmental factors, as reported for several softwood species [15,[31][32][33][34].
As well as annual ring width, the area-weighted basic density values significantly differed between sites (Table 3).In the present study, the quadratic function was well fitted to radial variation of basic density in relation to annual ring number from pith.The intercept parameter (β 2 ) for Kiso was higher than that for Hitachi (Table 5), implying that the mean basic density was higher at Kiso than that of Hitachi.In addition to the intercept parameter (β 2 ), the radial growth patterns of basic density are characterized by a combination of the shape (β 0 ) and variability (β 1 ) parameters (Fig. 2).
As the shape parameter (β 0 ) increases, the radial variation becomes more curvilinear.On the other hand, the smaller absolute value of variability parameter (β 1 ), the larger range of variability of radial variation.The higher shape parameter (β 0 ) was found for Hitachi, suggesting that basic density at Hitachi decreases more curvilinearly near the pith, followed by a minimum value and then a slight increase.The variability of basic density at Hitachi was slightly smaller than that of Kiso, because the lower variability parameters (β 1 ) were obtained for Hitachi.Based on both parameters (β 0 and β 1 ), it was concluded that the radial variations were more stable at Kiso than at Hitachi.In softwood species, wood density is often negatively affected by the growth rate [12].As shown in Table 5 and Fig. 3, radial growth rate at juvenile growth stage was faster at Hitachi than at Kiso.Differences between sites of basic density might be related to the differences in radial growth rate between sites.However, the slight difference between sites of basic density did not almost accompany difference in the boundary cambial age between juvenile and mature wood (Fig. 3).The boundary cambial age between juvenile and mature wood is considered to be around 13 annual ring number form pith, regardless of differences of the environment factors between sites.Therefore, the radial variation of basic density was somewhat affected by radial growth rates influenced by environmental factors, but the results suggest that there is little impact on final utilization as construction materials.
The best models for radial variation of basic density in both sites included the random effect of family but not of block (Table 4), suggesting that basic density is more affected by genetic factors than micro-environmental factors.The replicate blocks of two progeny test sites used in the present study were arranged from the upper to the lower part of the same slope.Therefore, genetic improvement of basic density is expected at slopes as steep as the progeny test sites investigated in the present study, without being affected by slope.In addition, the random effects in selected models were included in rate and intercept terms.This result implies that mean values and variation range from pith to bark of basic density were varied among families.Therefore, the mean basic density might improve by selecting the families with higher b 2 and lower b 1 , respectively.
To clarify the genetic control of radial growth in hinoki cypress, heritability and correlation between sites for annual ring width were estimated at each radial position (Table 6).The trend to increase heritability of annual ring width with age in the Hitachi site (Table 6) might be attributed to the widening with ages of the difference between families that acquire a better light environment at younger ages and families that did not.The low heritability of annual ring width at Kiso (Table 6) is thought to be due to the strong influence of environmental factors, including slope, soil type, snow, and animal damage.In addition, low correlations between sites in almost all radial positions except for 1st-3rd and 10th-12th radial positions (Table 6) might be caused by the difference in radial growth patterns between sites as shown in the results of modeling for radial growth.To maximize genetic gains and economic benefits from forest plantations, the strategy of selecting individuals for specific and known environments should be applied [14].Further study is needed to identify environmental factors that In several softwood species, wood density is generally known as a more heritable trait than radial growth [1,3,4,20].As well as previous reports on other softwood species, genetic control on basic density in hinoki cypress was stronger than that on radial growth.In addition, correlations between sites of basic density at all radial positions were significant and high.Moreover, the correlations between sites of random effects of family (b 1 and b 2 ) in the selected models were also estimated, resulting in significant correlations between sites was found in b 2 but not in b 1 .The low correlation in random effect b 2 between sites indicates that the range of basic density variation is affected by the interaction between genetics and environments.This result may be related to the difference between sites of radial growth patterns.In summary, although the range of variation may vary somewhat due to radial growth patterns, basic density in hinoki cypress can be effectively improved by tree breeding regardless of the environments in plantations such as climatic conditions.

Conclusion
To clarify the effects of genetics, environments, and interaction of those on radial growth and radial variation of basic density in hinoki cypress, we attempted to model the radial growth and radial variation of basic density using 36 hinoki cypress families planted at two progeny test sites with different climatic conditions.In comparing the fixed effect parameters of selected models for each site, the radial growth patterns slightly differed between sites.In addition, the stem diameter reaching the plateau might be varied among the progeny test's blocks.On the other hand, the radial variation patterns of basic density were slight differences between sites.Differences between sites of basic density might be related to the differences in radial growth rate between sites.However, radial variation of basic density is largely affected by genetic factors rather than micro-environmental factors.In addition, genetic control on basic density in hinoki cypress was stronger than that on radial growth, and the effects between genetic and environmental factors were smaller.In conclusion, although radial growth in hinoki cypress varies by the effects of micro-and macro-environmental factors and has some influences on radial variation of basic density, basic density is more strongly controlled by genetic factors, enabling effective improvement by tree breeding programs.

Fig. 3
Fig. 3 Cumulative annual ring width and basic density in relation to annual ring number from pith.CRW, cumulative annual ring width; BD, basic density.The circles in grey color indicate values in each sample tree.Regression curves indicate the formula with the fixed-effects parameters (Table5) and random effects of family or block in sites for the best-fitted models selected by AIC (Table4).The vertical solid lines indicate the annual ring number from the pith in which the rate of change in basic density reached less than 1% (13.9 years at Hitachi and 12.5 years at Kiso)

Fig. 4
Fig. 4 Correlations of random-effects parameter estimates of families in the best model for basic density between sites.r, correlation coefficient between sites.The values in parenthesis represent p-values

Table 1
Outline of two progeny test sitesOf six replicate blocks in the Kiso site, five replicate blocks were subjected for the present study because the number of survival individuals was limited in one replicate block.The climate dataset was based on a 1-km mesh scale via the digital national land information of Japan[21]

Table 2
Developed nonlinear models for cumulative annual ring width and basic density CRW ijk and BD ijk are cumulative annual ring width and basic density of the ith annual ring number from the pith of the jth individual in the kth family or block in progeny.RN ijk is the ith annual ring number of the jth individual of the kth family or block in progeny.β 0 , β 1 , and β 2 are the fixed-effect parameters.b 0jk , b 1jk , and b 2jk are the random effect of the jth individual nested by the kth family or block in progeny.e ijk is the residual

Table 3
Mean and standard deviation of measured traits and the result of t-test in measured traits between sites n, number of individuals; DBH, stem diameter at breast height; ARW, mean annual ring width of all annual ring width; BD, Area-weighted basic density; SD, standard deviation.The p-values were obtained by t-test between progeny test sites

Table 4
The AIC values of models for cumulative annual ring width and basic density CRW, cumulative annual ring width; BD, basic density.-, the model failed to converge or the model was not applied.Bold values indicate the minimum AIC values among developed models.Model forms refer to Table 2

Table 5
Fixed-effects parameters of the selected optimal model for radial growth pattern and basic density SE, standard error; CRW, cumulative annual ring width; BD, basic density.Model forms refer to Tables2and4

Table 6
Heritability and correlation between sites for annual ring width and basic density h 2 , narrow-sense heritability; r, the correlation coefficient between sites.The values in parentheses represent the p-value.The correlation coefficients in bold style indicate significance at a 5% level