Wood structure explained by complex spatial source-sink interactions

Overall framework

Cells in our model are arranged along independent radial files, with each cell in one of either the proliferation, enlargement-only, wall thickening, or mature zones, depending on the distance of the cell’s centre from the inside edge of the phloem and the time of year. Only cells that contribute to the formation of xylem tracheids are treated explicitly. A daily timestep is used, on which cells in the proliferation and enlargement-only zones can enlarge in the radial direction if these zones are non-dormant, and on which secondary-wall thickening can occur in the wall thickening zone. Cells in the proliferation zone divide periclinally if they reach a threshold radial length. Cell-size control at division is intermediate between a critical size and a critical increment²². Mother cells divide asymmetrically, with the subsequent relative growth rates of the daughters inversely proportional to their relative sizes. Size at division and asymmetry of division are computed with added statistical noise²², and therefore the model is run for an ensemble of independent radial files with perturbed initial conditions.

Equations and parameters

Cell enlargement and division

Cells in the proliferation and enlargement-only zones, when not dormant, enlarge in the radial direction at a rate dependent on temperature and relative sibling birth size. A Boltzmann-Arrhenius approach is used for the temperature dependence³⁰:

where μ is the relative rate of radial cell growth at temperature T (μm μm⁻¹ day⁻¹), μ₀ is μ at temperature T₀ (=283.15 K), E_a is the effective activation energy for cell enlargement, k is the Boltzmann constant (i.e. 8.617 x 10⁻⁵ eV K⁻¹), and T is temperature (K). μ₀ was calibrated to an observed mean radial file length at the end of the elongation period dataset²³ (Table 1; see “Observations”), and E_a was calibrated to an observed temperature dependence of annual ring width dataset³¹ (Table 1; Supplementary Fig. 4; see “Observations”).

Radial length of an individual cell then increases according to:

where ΔL_r is the radial increment of the cell (μm day⁻¹), L_r is the radial length of the cell (μm), and ϵ is the cell’s growth dependence on relative birth size, given by²²:

where g_asym is the strength of the dependence of relative growth rate on asymmetric division (Table 2; unitless), and α_b is the degree of asymmetry relative to the cell’s sister²² (scalar):

where L_r,_b is the radial length of the cell at birth (μm) and ({L}_{r}{,}_{b}^{sis}) is the radial length of its sister at birth (μm), which are calculated stochastically²²:

where L_r,_d is the length of the mother cell when it divides (μm) and Z_a is Gaussian noise with mean zero and standard deviation σ_a (Table 2; −0.49 ≤Z_a≤ 0.49 for numerical stability).

Length at division is derived as²²:

where f is the mode of cell-size regulation (Table 2; unitless), χ_b is the mean cell birth size (Table 3; μm), and Z is Gaussian noise with mean zero and standard deviation σ (Table 2).

The first cell in each radial file is an initial, which produces phloem mother cells outwards and xylem mother cells inwards. It grows and divides as other cells in the proliferation zone, but on division one of the daughters is stochastically assigned to phloem or xylem, the other remaining as the initial. The probability of the daughter being on the phloem side is f_phloem (Table 3).

Cell-wall growth

Both primary and secondary cell-wall growth are influenced by temperature, carbohydrate concentration, and lumen volume. A Michaelis-Menten equation is used to relate the rate of wall growth to the concentration of carbohydrates in the cytoplasm:

where ΔM is the rate of cell-wall growth (mg cell⁻¹ day⁻¹), ΔM_max is the carbohydrate-saturated rate of wall growth (mg cell⁻¹ day⁻¹), θ is the concentration of carbohydrates in the cell’s cytoplasm (mg ml⁻¹), and K_m is the effective Michaelis constant (mg ml⁻¹; Table 1).

The maximum rate of cell-wall growth, ΔM_max, is assumed to depend linearly on lumen volume (a proxy for the amount of machinery for wall growth), and on temperature as in Eq. (1):

where ω is the normalised rate of cell-wall mass growth (i.e. the rate at T₀; Table 1; mg ml⁻¹ day⁻¹), V_l is the cell lumen volume (ml cell⁻¹), and E_aw is the effective activation energy for wall building (eV; Table 1). ω and K_m were calibrated to an observed distribution of carbohydrates²³ (see next section). E_aw was calibrated to an observed temperature dependence of maximum density³¹ (Table 1; see “Observations”).

Lumen volume is given by:

where V_c is total cell volume (ml cell⁻¹) and V_w is total wall volume (ml cell⁻¹). Cells are assumed cuboid and therefore V_c is given by:

where L_a is axial length (μm; Table 2) and L_t is tangential length (μm; Table 3). V_w is given by:

where M is wall mass (mg cell⁻¹) and ρ is wall-mass density (Table 2; mg[DM] ml⁻¹).

Cells in the proliferation and enlargement-only zones only have primary cell walls. ΔM_max (Eq. (9)) is therefore given the following limit:

where ({V}_{{w}_{p}}) is the required primary wall volume:

where W_p is primary cell-wall thickness (Table 3; μm).

Carbohydrate distribution

The distribution of carbohydrates across each radial file is calculated independently from the balance of diffusion from the phloem and the uptake into primary and secondary cell walls. The carbohydrate concentration in the phloem is prescribed at the mean value observed across the three observational dates in²³, as described below in “Observations”, and the resulting concentration in the cytoplasm of the furthest living cell from the phloem is solved numerically. The inside wall of this cell is assumed to be impermeable to carbohydrates and therefore provides the inner boundary to the solution. It is assumed that the rate of diffusion across each file is rapid relative to the rate of cell-wall building, and therefore concentrations are assumed to be in equilibrium on each day. Carbohydrate diffusion between living cells is assumed to be proportional to the concentration gradient:

where q_i is the rate of carbohydrate diffusion from cell i − 1 to cell i (mg day⁻¹) and η is the resistance to flow between cells (calibrated to the observed distribution of carbohydrates²³, see next section; Table 1; day ml⁻¹).

As it is assumed that carbohydrates cannot diffuse between radial files, at equilibrium the flux into a given cell must equal the rate of wall growth in that cell plus the wall growth in all cells further along the radial file away from the phloem. From this it can be shown that the equilibrium carbohydrate concentration in the furthest living cell from the phloem in each radial file is given by:

where θ_p is the concentration of carbohydrates in the phloem (Table 3; mg ml⁻¹) and n is the number of living cells in the file (phloem mother cells are ignored for simplicity). The rate of wall growth in each cell depends on the concentration of carbohydrates (Eq. (8)), and therefore θ_n must be found that results in an equilibrium flow across the radial file. This is done using Brent’s method⁴¹ as implemented in the “ZBRENT” function⁴².

Zone widths

The widths of the proliferation, enlargement-only, and secondary wall thickening zones vary through the year, and are fitted to observations on three dates²³ (see Supplementary Fig. 2 and “Observations”). Linear responses to daylength were found, which are therefore used to determine widths for the observational period and later days:

where z_k is the distance of the inner edge of the zone from the inner edge of the phloem (μm), k is proliferation (p), secondary wall thickening (t), or enlargement-only (e), a_k and b_k are constants (Table 3), dl is daylength (s), and DOY is day-of-year. The proliferation zone width on earlier days when non-dormant was fixed at its DOY 185 width (assuming this to be its maximum, and that it would reach its maximum very soon after cambial dormancy is broken in the spring). During dormancy, the proliferation zone width is fixed at its value on DOY 231 (the first day of dormancy²³). The enlargement-only zone width prior to DOY 185, the first observational day, is assumed to be a linear extension of the rate of change after DOY 185. The wall-thickening zone width plays little role prior to DOY 185 at the focal site, and so was set to its Eq. (17) value each earlier day. On all days the condition z_t≥z_e≥z_p is imposed, and zone widths cannot exceed their values at 24 h daylength (necessary for sites north of the Arctic circle). Supplementary Figure 2 shows the resulting progression of zone widths through the year, together with the observed values.

Dormancy

Proliferation was observed to be finished by DOY 231²³, and so the proliferation and enlargement-only zones are assumed to enter dormancy then. Release from dormancy in the spring is calculated using an empirical thermal time/chilling model³³. It was necessary to adjust the asymptote and temperature threshold of the published model because the heat sum on the day of release calculated from observations in Sweden (see “Observations”) was much lower than reported for Sitka spruce buds in Britain in the original work:

where dd_req is the required sum of degree-days (°C) from DOY 32 for dormancy release and cd is the chill-day sum from DOY 306. The degree-day sum is the sum of daily mean temperatures above 0 °C, and the chill-day sum is the number of days with mean temperatures below 0 °C. Dormancy can only be released during the first half of the year.

Simulation protocols

Each simulation consisted of an ensemble of 100 independent radial files. Each radial file was initialised by producing a file of 100 cells with radial lengths χ_b(1+Z_a), allowing these to divide once, ignoring the second daughter from each division, and then limiting the remaining daughters to only those falling inside the proliferation zone on DOY 1. Values for ϵ (the relative growth of daughter cells) and L_r,_d (the radial length at division) were derived for each cell. The main simulations were conducted at the observation site in boreal Sweden (64.35°N, 19.77°E) over 1951–1995 to capture the growth period of the observed trees²³. Results are mostly presented for 1995 when the observations were made. Simulations for calibration of the effective activation energies (i.e. E_a and E_aw) were performed at 68.26°N, 19.63°E in Arctic Sweden over 1901–2004³¹. Daily mean temperatures for both sites were derived from the appropriate gridbox in a 6 h 1/2 degree global-gridded dataset⁴³.

Observations

Observations of cellular characteristics and carbohydrate concentrations²³ were used to derive a number of model parameters, and to test model output (model calibration and testing were performed using different outputs). According to the published study we used, samples were cut from three 44 yr old Scots pine trees growing in Sweden (64°21’ N; 19°46’ E) at different times during the growing season. 30 μm thick longitudinal tangential sections of the cambial region were made, and the radial distributions of soluble carbohydrates measured using microanalytical techniques²³. Cell sizes, wall thicknesses, and positions in their Fig. 1²³, an image of transverse sections on three sampling dates, were digitised using “WebPlotDigitizer”⁴⁴. These, together with the numbers of cells in each zone and their sizes given in the text of that paper, were used to estimate zone widths, which were then regressed against daylength to give the parameters for Eq. (17) (Table 3), mean cell size in the proliferation zone on the first sampling date (used to derive χ_b; Table 3), mean cell tangential length (Table 3), and final ring width (used to calibrate μ₀; Table 1). The thickness of the primary cell wall (Table 3) was derived by plotting cell-wall thickness against time and taking the low asymptote.

The distributions of carbohydrates along the radial files on the last sampling date for “Tree 1” and “Tree 3” (results for “Tree 2” were not shown for this date) shown in Fig. 2 of the observational paper²³ were calculated. The masses for each of sucrose, glucose, and fructose in each 30 μm section were digitised using the same method as for cell properties and then summed and converted to concentrations, with the results shown in Supplementary Figure 5. Mean observed carbohydrate concentrations and cell masses at four points were used to calibrate values for the η, ω, and K_m parameters in Table 1. Calibration was performed by minimising the summed relative error across the observations.

The calibration target for the effective activation energy for wall deposition (i.e. E_aw) was the observed relationship between maximum density and mean June-July-August temperature over 1901-2004 in northern Sweden³¹ (Supplementary Fig. 3), and for the effective activation energy for cell enlargement the relationship between ring width and temperature (i.e. E_a) target was the same study (Supplementary Fig. 4). These observations were made on living and subfossil Scots pine sample material from the Lake Tornesträsk area (68.21–68.31°N; 19.45–19.80°E; 350–450 m a.s.l.) using X-ray densitometry for maximum density, and standardised to remove non-climatic information³¹.

Reporting summary

Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.

Source: Ecology - nature.com