^{2} values increasing from around 0.6 to 0.7 in the model validation. Bias and RMSE were also lower for the generalized function. However, error analysis showed that there were no large differences between these two methodologies, evidencing that GWR provides results which are as good as the more frequently used LMM methodology, at least when no additional measurements are available for calibrating.

In forestry, spatial heterogeneity is theorized as one of the major drivers of biological diversity (

Global and local regression models have traditionally been used for studying spatially heterogeneous data (

There are several methods for solving the heterogeneity problem, for example, models based on heteroscedasticity specification, generalized additive models, classification and regression trees, linear mixed models (LMM) and geographically weighted regression (GWR). In this document the latter two methods were analysed.

LMM has been widely used for height-diameter modelling (

The second methodology analysed, GWR, is used to explore spatially heterogeneous processes (

Tree height is an important variable for forest management in several ways: (1) total height together with diameter at breast height allows us to estimate the total volume of the tree, (2) in regular stands, dominant height together with age are used to estimate the site index and (3) the distribution of tree heights in the stand is very important to describe the vertical structure.

Nevertheless, height measurement can be time consuming and consequently expensive, therefore it is common to use height-diameter models instead. There are two kinds of height-diameter models: local models and generalized models.

Local models are designed specifically for forest stands and are the most precise curves for a particular forest stand (

The main aim of this paper is to analyse the suitability of GWR for modelling height-diameter relationships in a large area and to compare it with LMM, which is one of the most frequently used methods for including spatial heterogeneity in these models.

To test whether GWR provides as good results as the more frequently used LMM, we applied both methodologies to two different functions: the

To study GWR, we selected mono-specific stands of Aleppo pine (

Data were obtained from 230 sample plots with 2582 trees belonging to the Third Spanish National Forest Inventory (SNFI). SNFI design is based on permanent sample plots located at the nodes of one kilometre square grid, re-measured in an inventory cycle of ten years. The sample plots consisted of four concentric circles of 5, 10, 15 and 25 meter radius, in which the diameters and heights of all trees with a breast height diameter over 7.5, 12.5, 22.5 and 42.5 cm respectively were measured. Furthermore, plot centre coordinates, together with the polar coordinates of each measured tree, i.e. distance and azimuth form the plot centre, were recorded (

Data in the studied area were randomly split into two subsets. The first consisted of 152 plots with 1697 trees and was used to develop the model (training data). The second was used to validate the model (validation data) and consisted of 78 plots with 885 trees (

Two different height-diameter functions were selected for the study. The first of these was the local height-diameter Petterson’s function (

where _{0} and_{1} are the parameters to be estimated.

The second selected function was the Cañadas I generalized height-diameter function (

where _{o} is dominant height (m) obtained according to the Hart criteria (_{0} is dominant diameter, calculated as the average of the 100 thickest trees per hectare, and _{0} is the parameter to be estimated.

As GWR is basically a technique developed for linear regressions (

Petterson:

Cañadas I:

These two height-diameter functions were fitted following the two different methodologies detailed below.

In a first step, linear mixed models were used for fitting the height-diameter equations. As data came from a hierarchical structure, with several trees measured in the same plot, a mixed model was used with the plot as the grouping structure of the random effects. We included random effects in the _{0}_{ }and _{1}_{ }term of Petterson local model (_{0}_{ }for the Cañadas I (

The models were fitted using the maximum likelihood method of the lme procedure (

In a second step, the models were fitted using the GWR with three different bandwidths and an optimal number of trees.

Two of these bandwidths, calculated taking into account the SNFI grid, were 1000 m and 500 m. These two values were used to determine whether the SNFI design influenced the GWR results. The third bandwidth was calculated according to the calibration recommended by

Gaussian weighted function:

where _{ik} is the _{ik} is the distance between observations

Finally, the optimal number of trees was calculated according to the calibration recommended by

Bisquare weighted function:

where _{ij} is the _{ij} is the distance between observations

GWR provides a range of coefficients geographically located in the study area, and a continuous GWR map of coefficients can be obtained by the coefficients interpolation. The Inverse Distance Weighting (

GWR was fitted using the GWModel package in R language (

Models were validated by applying GWR as well as LMM to the validation dataset, obtaining a set of predicted heights. In the case of GWR, the predicted values were obtained through the GWR map of coefficients.

In order to validate LMM, it was necessary to calibrate the model. Very few measurements can be used in the calibration, even less than the number of random parameters to be predicted (

The validation was performed by comparing observed and predicted heights through a statistical analysis of the residuals. Mean error (ME), mean absolute error (MAE), root mean square error (RMSE) and coefficient of determination (R^{2}) were calculated according to the formula in

^{2}) and the Akaike information criterion (AIC).

For the Petterson height-diameter function, the coefficient _{0} mean value varied between 0.37 and 0.55 for LMM and also for GWR independently of the bandwidth (between 0.34 and 0.52 for a bandwidth of 1000, between 0.32 and 0.69 for a bandwidth of 500, between 0.34 and 0.71 for a Fixed Kernel bandwidth, and between -0.19 and 0.64 for an adaptive kernel bandwidth). The range of this parameter was lowest when obtained by LMM or by GWR with a bandwidth of 1000 m, but increased for the bandwidth of 500 and the fixed kernel, and reached a maximum for the adaptive kernel. The coefficient a_{1,} i.e. the slope of the curve, showed higher variation with larger standard deviations and wider ranges than the a_{0}. So, for LMM the standard deviation of a_{1} was 0.48 while for GWR were from 0.65 to 1.39. However, mean values along the study area were similar for all the methodologies, from 1.79 to 2.06. Once again, the largest variations were found for the GWR fitted with fixed or adaptive kernel, proving to be more sensitive to the spatial variability of the data. Moreover, the smallest AIC was obtained with the fixed kernel, representing an optimal bandwidth of 363 m; and with the adaptive kernel, representing an optimal number of trees of 26. However, all GWR methods showed certain incongruent results such as negative or very large parameter values. An example of estimation using the Petterson map of coefficients is presented in

In the case of the Cañadas-I height-diameter function it can be seen that the coefficient _{0} had slightly larger standard deviations and ranges for GWR than for LMM. Additionally, the mean values were almost the same between both methods. The AIC for this function was smaller when GWR was fitted with a predefined bandwidth of 500 m. For this generalized height-diameter function, the optimal bandwidth calculated by fixed kernel was 755 m, and the optimal number of trees in the adaptive kernel was 65. Again, some incongruent minimum values can be seen for the GWR methodology in

As response variables in the Petterson and Cañadas-I functions were different, the AIC and R^{2} values were presented in

The validation revealed that the Petterson height-diameter function had poorer predictive capacity than the Cañadas I generalized function. Consequently, it can be seen that the R^{2 }for the LMM with only fixed effects increased from 0.27 to 0.71, when calibration the improvement was from 0.63 to 0.72 and for GWR the R^{2 }ranged between 0.29 and 0.35 for Petterson while for Cañadas was 0.71.

When analysing the errors obtained from Petterson model, it was observed that, independently of bandwidth or kernel, GWR resulted in lower bias than LMM with only fixed effects (ME=0.6) or in similar than the calibrated (ME=0.1). Moreover, MAE and RMSE were very similar for GWR and LMM with only fixed effects (MAE=1.3 and RMSE=1.7 for LMM without calibration; and MAE=1.3 - 1.4 and RMSE=1.7 for GWR). However both MAE and RMSE were smaller when calibrating LMM (MAE= 0.9 and RMSE=1.3). Regarding the Cañadas-I height-diameter function, the error analysis showed similar results for all the GWR and LMM methods with and without calibration (

The standardized residual plots showed that residuals are symmetrically distributed in a random pattern around the horizontal axis for any of the functions and fitting methodologies (

Tree height is an important variable which is used for estimating stand volume or site quality (

In the case of the Petterson model, it was observed that this method had insufficient predictive quality when GWR or uncalibrated LMM was applied. It is possible that the large heterogeneity of the studied area made it difficult for the Petterson function to reach the necessary quality. This could be due to the fact that the Petterson formula only had diameter at breast height as an independent variable and would appear to be too simple to reproduce the forest heterogeneity. Nevertheless, when Petterson equation was validated with calibrated LMM, a high improvement in the predictive quality was observed. In the case of height-diameter model, random component represents the specific conditions of every stand in the model. It causes better precision of the model for particular location in comparison with the model with fixed parameters. This could be important in a local model such as Petterson, designed specifically for giving precise prediction in individual forest stands (

In contrast, the Cañadas I model revealed high stability in the model fitting and high predictive quality in the validation in both LMM and GWR. In these cases, R^{2} values in the validation improve from local to generalized regressions from 0.63 to 0.72 in calibrated LMM. Similarly, MAE and RMSE were lower for the generalized models. The comparison between local and generalized equations is analogous with previous studies which concluded that the inclusion of stand-density variables such as dominant height or dominant diameter in the base height–diameter function increased the accuracy of prediction (

The predictive quality of the model was improved when few additional height measurements were used for estimating the random parameters for the validation dataset. Calibration using small subsamples of trees is one of the main utilities of mixed models in forest applications (

In some of the GWR models the presence of odd, extreme values was observed. With the Petterson equation these values appear in every GWR model and did not correspond to specific plots but rather, were scattered randomly in the study area. This could be another indicator of the weakness of these models. However, in the case of Cañadas I, these odd values only appeared in the models with a bandwidth of 500 m and adaptive kernel. Specifically, this odd value was located in plot number 1473, which contained only 3 trees and was located in a relatively isolated area. Odd values have also been observed in similar height-diameter modelling studies (

In our study, the observed GWR models gave the same residuals regardless of the calibration used. However, according to

The maps of coefficient in _{0} takes the highest values and a_{1} the lowest; therefore in

When comparing the LMM and GWR models in this study, it was found that GWR estimation produces similar errors to the LMM when generalized equation was applied. Analogous results have been found in other studies, supporting our finding that GWR gives similar quality estimations to the LMM methodology and predicts the response variable accurately (

LMM and GWR provide different perspectives and insights for data analysis and model prediction. In LMM, the fixed component explains the impact of different variables as with the ordinary least squares regression (

Two different methodologies were used to model local and generalized height-diameter relationships: the geostatistical method GWR and LMM methodology. Additionally, two different height-diameter functions were selected to fit these models: a local function and a generalized function. Finally, the quality of the estimated models was compared through statistical analysis.

As expected, it was observed that in general, both LMM and GWR provide better prediction quality when applied to a generalized height-diameter equation. However, error analysis revealed that there were no great differences between these two methodologies, evidencing that GWR applied to a generalized function provides results which are at least as good as those obtained using the calibrated LMM methodology. However, GWR has a number of advantages over LMM. For instance, GWR provides information on a continuous map of regression coefficients within the studied area. Consequently, is possible to calculate tree heights in every single location without the need for additional measurements. Furthermore, the GWR method can help us to understand spatial processes and may have interesting applications in the area of forest yield and production.

The authors would like to thank two anonymous reviewers and handling editor for their valuable comments that helped improve the quality of this study.