Graft survival of Pinus engelmannii Carr. in relation to two grafting techniques with dormant and sprouting buds

Rootstock production

The rootstocks were produced in the “General Francisco Villa” forest nursery in Durango, Mexico (23° 58′ 20″ North and 104° 35′ 56″ West, at 1,875 m elevation). During the first year, the rootstocks were grown in polystyrene trays with 77 cavities (each with a capacity of 170 mL); the substrate was a mixture of equal parts of peat moss and composted pine bark. In the second year, each rootstock was transplanted into a 3.5 L black polyethylene bag, which contained equal parts of forest soil and pine bark.

The main physical properties of the substrate were as follows: electrical conductivity, 0.16 mmhos cm⁻¹; pH 6.04; organic matter content, 3.8%; and organic carbon content, 2.2%. The total porosity, aeration porosity and the water-holding capacity of the substrate were calculated with the methodology proposed by Landis et al. (1990), yielding values of 82%, 44% and 38% respectively. Analysis of the chemical properties of the substrate revealed values of 1.02% for nitrogen-ammonium (N-NH4), 16.73% for phosphorus (P) and 1.8% for potassium (K).

At age four years (July 2018), and six months before grafting, the rootstocks were transplanted to 5 L containers (bags), to favour reactivation of the root system. The substrate was forest soil, and to each container we added 50 g of Multicote^TM 6: 18-6-12 (N-P-K) + micronutrients, which is a slow-release fertilizer manufactured by the Haifa Group^TM, with a longevity of six months at an average temperature of 21 °C. The rootstocks were watered every three days, with two L of water per plant.

Collection of scions

The scions were taken from trees in a commercial plantation of P. engelmannii, in the Ejido Aquiles Serdán, Durango, Mexico (location 23° 53′ 39″ North and 104° 33′ 44″ West and 1,898 m elevation). The plantation was five years old when the scions were collected and showed good adaptation and growth. The selected donor trees were of average height three m and diameter 7 cm at the base of the stem. The scions were collected the day before grafting and placed in layers in 72 L plastic boxes; each layer of scions was covered with sawdust wetted with a solution of three g L⁻¹ of Captán^® fungicide, to prevent fungal damage.

Grafting

Grafting was carried out in the nursery of the Institute of Forestry and Wood Industry of the Universidad Juárez del Estado de Durango (ISIMA-UJED), in a greenhouse of dimension 6 × 8 × 3 m (width × length × height), covered with white plastic (caliber 720). To prevent excessive temperatures being reached and to promote adequate relative humidity, two shading meshes (providing 50% and 70% cover) were placed over the greenhouse, 30 cm above the plastic cover. In addition, two air conditioning systems were installed at each end of the greenhouse. The average temperature was maintained at 22 °C (maximum, 26 °C and minimum, 7 °C) and the relative humidity fluctuated between 72% and 82%.

Treatments evaluated

The grafting methods evaluated were the side veneer technique described by Muñoz et al. (2013) and Pérez-Luna et al. (2019) and the top cleft technique described by Muñoz et al. (2013). For each technique, two phenological stages of the buds were tested; half of the grafts were made with scion buds at the end of latency (grafted on January 18, 2019) and the other half were made with the scion buds at the beginning of the sprouting stage (grafted on February 21, 2019). In total, four treatments were evaluated, with 120 grafts (60 of side veneer and 60 of top cleft). Thirty grafts of each phenological condition of the buds were made with each grafting technique.

For the side veneer graft, a longitudinal cut of six cm was made on one side of the scion and a wedge of approximately one centimetre was formed with a cut at the lower end of the other side (Fig. 1A). A lateral cut was made in the rootstock, of the same length as the cut on the scion, and a one-cm slit was left at the end of the cut on the rootstock, into which the scion wedge was inserted and tied (Fig. 1B).

For the top cleft graft, two longitudinal 6 cm cuts were made at the bottom and opposite sides of the scion, forming a wedge shape; the central leader was eliminated from the rootstock, and a central cut (fissure) was made to an approximate depth of 6 cm, into which the scion wedge was inserted (Fig. 1C). The two components of the graft were then tied together (Fig. 1D).

The grafts were tied (both techniques) with transparent rubber tape and sealed with vinyl paint mixed with 3 g L⁻¹ of Captán^® fungicide, to prevent pathogens entering the graft union. Finally, a 5 L transparent plastic bag was placed around the grafted area, into which 1 L of water was poured to generate a humid microenvironment (Fig. 1E); in addition, each graft was covered with a kraft paper bag, which provided protection from solar radiation (Fig. 1F).

The grafts were watered every three days with plain water. In addition, from the 3rd month of evaluation, Promyl^® fungicide was added to the irrigation water (two g L⁻¹) every eight days, to prevent fungal damage. To compensate for possible nutrient deficiency due to the loss of mycorrhiza caused by fungicide application, fertigation was applied during the six-month evaluation period by adding Ultrasol^® Triple 19 (N-P-K) + MgO + micronutrients, a water-soluble fertilizer manufactured by Soquimich commercial, at a dose of 3 g L⁻¹ to the water. The solution was applied in a 7 L manual watering can.

Experimental design and variables evaluated

The treatments were applied in a 2 × 2 factorial experimental design (two grafting techniques × two phenological stages of the scion buds). Each experimental unit consisted of 10 grafts, with three repetitions per treatment. Before grafting, seven internal grafting variables of the scion and rootstock were measured (Table 1).

Sprouting and survival of the grafts were evaluated monthly for six months, and the data were examined by the Kolmogorov–Smirnov test to evaluate the normality of the variables evaluated. Analysis of variance (ANOVA) was used to detect potential significant differences between treatments. When significant differences were indicated and the variables were normally distributed, a post hoc Tukey’s means test was carried out, with an initial confidence interval of 95% (α = 0.05). The Bonferroni correction was applied to reduce the probability of making a type I error (Garamszegi, 2006; Napierala, 2012); the corrected significance value was α* = 0.0125. In addition, Student’s t-tests were used to determine any significant differences in graft survival due to different levels of the independent factors, i.e., side veneer grafts vs top cleft grafts and scions with buds at the end of latency vs scions with buds at the beginning of sprouting. In order to reduce the effect of extreme observations, before performing the analysis of variance and Student’s t tests (Burbidge, Magee & Robb, 1988) the survival value of each treatment was transformed, by calculating the square root of the sine function of the survival quotients.

Finally, the Weibull’s accelerated failure time model was fitted to the data to predict the estimated graft survival time after the evaluation period, and the associated hazard ratio was also calculated. Other hazard ratios were estimated using the Cox’s proportional hazards model. All statistical analyses were implemented in the “Survival” package in the free R software (R Development Core Team, 2018).

Fitting the Weibull’s accelerated failure time model and the hazard ratio

To fit the accelerated failure time model to the graft survival data, the following independent variables were used: grafting technique, phenological stage of the scion buds and the seven grafting variables (Table 1). To detect the significant variables affecting the estimated survival time of the grafts, stepwise regression was applied using the “StepReg” package in R (R Development Core Team, 2018). The Weibull´s accelerated failure time model is defined as follows:

(1) $$Ln\left(T\right)=\mathrm{\alpha }+\delta x_i+\sigma \epsilon$$

where Ln (T) is the natural logarithm of the mean survival time (T) after the study, i.e., at time (T) at least one death may occur among the grafts that were alive at the end of the evaluation period; α is a scale parameter of the model, δ is the coefficient of the explanatory variable, σ is the shape parameter of the model and ε is the error of the distribution function (George, Seals & Aban, 2014; Zhang, 2016). The Weibull’s hazard ratio is described as follows:

(2) $$H{R}_W={\left(e^{-\mathrm{\text{β}}}\right)}^{\mathrm{\lambda }-1}$$

where HR_W is the Weibull’s hazard ratio, which represents the increase or decrease in the risk of death of grafts as a function of the independent variables; β is a coefficient that indicates the effect of a given independent variable (from the same xi variables used in Eq. (1)). The value of the dependent variable is not used in Eq. (2), because the Weibull hazard ratio is considered to be “constant for each variable”, i.e., there is only one hazard ratio value for each explanatory variable in the model (Igl, 2018). The β coefficient is interpreted as follows: if β is positive the risk of death decreases as the value of the independent variable decreases, and if β is negative, the risk of death decreases as the value of the independent variable increases (Carroll, 2003). Finally, λ is a shape parameter of the model. If λ > 1 the hazard ratio increases, and if λ < 1 the hazard ratio decreases (Zhang, 2016). The β and λ parameters were calculated with flexible parametric regression, which is useful for modelling survival from a time of origin to the instant at which an event occurs (life or death) (Igl, 2018; Pérez-Luna et al., 2020b).

The Weibull hazard ratio can take values ranging between 0 and ∞ (Ruíz, 2012). A hazard ratio of 0 indicates that the risk of death decreases by 100% [(0 − 1) × 100]. A hazard ratio of > 0 but < 1, for example HR_W = 0.5, indicates that the risk of death due to the effect of a variable is reduced by 50% [(0.5 − 1) × 100] when the level of that particular variable changes. On the other hand, if the hazard ratio is 1.0, the risk of death does not vary due to the effect of changes in some particular variable [(1 − 1) * 100]; if the hazard ratio takes a value of two the risk of death doubles, i.e., the risk of death increases by 100% due to the effect of changes in the explanatory variable [(2 − 1) × 100] (Hilsenbeck et al., 1998; Spruance et al., 2004). Interpretation of the hazard ratio is similar for values greater than two (Ruíz, 2012). It is important to take into account that the Weibull hazard ratio is calculated individually for each dependent variable in the model (Zhang, 2016).

For fitting these models, dummy variables were used to code the dependent variable (survival) and the independent variables (grafting technique and phenological stage). Therefore, a live graft was coded as zero and a dead graft as one (censor variables). The coding for the grafting technique was one for side veneer grafts and two for top cleft grafts; the phenological stage was coded as one for cuttings with buds at the end of dormancy and two for cuttings with buds at the beginning of sprouting. For more details on the Weibull’s accelerated failure time model and its hazard ratio, see Pérez-Luna et al. (2020b).

Fitting the Cox’s proportional hazards model

The variables shown in Table 1 were used to fit this model, and the most significant variables were selected by stepwise regression. The Cox’s proportional hazards model used for calculating hazard ratios is as follows:

(3) $$H{R}_C=e^{\left({\mathrm{\phi }}_1{x}_{i1}+\dots +{\mathrm{\phi }}_{\mathrm{k}}{x}_{ik}\right)}$$

where HR_C is the Cox’s hazard ratio expressed as a function of the independent variables (x_i) with which the model is fitted in order to predict the increase or decrease in the risk of death of the individuals under study (grafts). φ represents the fit parameters of the model, up to the k-th independent variable. If a given φ is positive, the hazard ratio increases when the value of the corresponding x increases (decreasing the probability of survival); if the value of a φ is negative, the hazard ratio decreases when the respective x increases (increasing the probability of survival). The range of values and the interpretation of the Cox’s hazard ratio are the same as for the Weibull’s hazard ratio. The Cox’s hazard ratio value is calculated globally, i.e., by including all the independent variables at one time in the $H{R}_C$ equation (Ata & Sözer, 2007).

Several authors recommend the use of hazard ratios derived from the Weibull’s model, as long as the shape parameter of the model (λ), which must be calculated before using the hazard ratio (Lee & Wang, 2003; Lawless, 2011), is known; however, its counterpart, the hazard ratio of the Cox’s proportional hazards model, does not depend on the evaluation time or predicted survival and therefore has the advantage of being less restrictive than the HR_W (Cox, 1972).

To calculate the hazard ratio for graft survival using the Cox’s proportional hazards model, dummy variables, including censor variables (dead grafts), were also used, so that the coding for the dependent variables was the same as in the Weibull model. Further details on the use of this model to assess graft survival are given by Pérez-Luna et al. (2019).

Survival

A total of six months after grafting, the average survival of the grafts made with the top cleft and side veneer techniques were 56.7% and 18.3%, respectively, with significant differences between them. There were also differences in regard to the phenological stages of the buds, analyzed as individual factors (Table 2 and Fig. 2), yielding 50% survival in grafts with buds at the end of latency and 25% survival in the grafts with buds at the beginning of sprouting.

The analysis of variance indicated significant differences (before the Bonferroni correction) due to the effect of these two treatments (Table 2); however, after Bonferroni correction only the grafting technique was statistically significant. Furthermore, the Tukey test revealed that the best interaction was the combination of top cleft grafts with the scion buds at the end of the latency, with 80% survival, while the treatment yielding lowest survival (16.7%) was the side veneer grafts with the scion buds at the beginning of sprouting (Fig. 3).

Fitting the Weibull’s accelerated failure time model and the hazard ratio

The stepwise regression indicated that survival can be described with only four independent variables, using the Weibull’s accelerated failure time model for which the lowest value of the Akaike information criterion (AIC) was achieved (Table 3). The model fit was highly significant, even after Bonferroni correction (p < 0.0001), for predicting the mean graft survival time after finishing the six-month evaluation period.

In addition, all parameters estimated for the fitting of this four-variable model were significant for the original critical level proposed at the beginning (p < 0.05), and only the phenological stage of the buds was not significant after the Bonferroni correction (p < 0.0125) (Table 4).

Application of the accelerated failure time model (Eq. (1)) to the data showed that the estimated time of survival, after finishing the evaluation period, was greater for those grafts produced by the top cleft technique (x = 2) and with scions’ buds at the end of latency (x = 1). In addition, the estimated survival time was greater when the diameter of the scion was less than 11.4 mm (x ≤ 1) and the height of the rootstock greater than 58.5 cm (x ≥ 1). The longest estimated average survival time, after the observation period (6 months), was 457 days, while the shortest estimated average survival time was 164 days.

The estimates depend on the combinations of the values of each independent variable, and these results represent the estimated time during which at least one death could occur among the grafts that were alive at the end of the evaluation period (Zhang, 2016; Pérez-Luna et al., 2020b).

The estimated parameters (calculated by regression) for fitting the Weibull hazard function (Zhang, 2016) are shown in Table 5. A negative sign of the β coefficient indicates that when the variable “grafting technique” increases by one unit, from x = 1 (side veneer graft) to x = 2 (top cleft graft), the risk of death of at least one graft will be reduced by 76% [(HR_W − 1) × 100 = (0.24 − 1) × 100 = −76%]; that is, the risk of death decreases by 76% when grafting with the top cleft technique. Interpretation of the negative value of the β estimator of the variable “rootstock height” is the same as the case described above; therefore, for rootstocks of height greater than 58.5 cm, the risk of death of the grafts is reduced by 4% [(HR_W − 1) × 100 = (0.96 − 1) × 100 = −4%]. On the other hand, the positive signs of β for the variables “phenological stage of buds” and “diameter of scion” indicate that by increasing the respective variable by one unit, the risk of death of the grafts will increase according to the value obtained for each hazard ratio; thus, in the case of the “phenological stage of buds” for x = 2 (grafts with buds at the beginning of sprouting), the risk of death increased by 64% [(HR_W − 1) × 100 = (1.64 − 1) × 100 = 64%] when the phenological stage is x = 1 (grafts with buds of scions at the end of dormancy). Finally, for scions of diameter greater than 13.38 mm, the risk of death of the grafts increased by 25% [(HR_W − 1) × 100 = (1.25 − 1) × 100 = 25%].

Fitting the Cox’s proportional hazards model

Using stepwise regression, the variables that best explained the hazard ratio in graft survival were also selected by fitting the Cox proportional hazards model and it was confirmed that the grafting technique, the phenological stage of the buds, the rootstock height and the scion diameter were the outstanding variables for the best fit (Table 6).

The variables that best explained the hazard ratio of the Cox’s proportional hazards model were the same variables selected for fitting the Weibull’s accelerated failure time model, and although Tables 3 and 6 are similar, note that the value of the Akaike’s information criterion was lower in Table 6 (AIC = 609.6) than that shown in Table 3 (AIC = 904.8). Therefore, the Cox’s proportional hazards model performed better than Weibull’s model for predicting the survival probability. The estimators calculated for fitting the Cox’s model are shown in Table 7.

The hazard ratio values (Table 8) were calculated using Eq. (3). It was estimated that the lowest value of the hazard ratio of grafts death (0.18) would be obtained for the top cleft grafts with the scion buds at the end of the dormancy period, using rootstocks taller than 58.5 cm and scions with diameter smaller than 11.4 mm; this value implies that the risk of death would decrease by 82% [(HR_C − 1) × 100 = (0.18 − 1) × 100 = −82%] when grafts are produced with this combination of variables, relative to the hazard ratio obtained for other combinations. On the other hand, the highest estimated hazard ratio corresponded to the side veneer grafts using scions with buds at the beginning of sprouting, for rootstocks shorter than 58.5 cm and scions with diameter greater than 11.4 mm, taking a value of 2.0, which indicates that the risk of death under this combination increases by 100% [(HR_C − 1) × 100; (2 − 1) × 100 = 100%], relative to the other combinations.