Title: | Flexible Tree Taper Curves Based on Semiparametric Mixed Models |
---|---|
Description: | Implementation of functions for fitting taper curves (a semiparametric linear mixed effects taper model) to diameter measurements along stems. Further functions are provided to estimate the uncertainty around the predicted curves, to calculate timber volume (also by sections) and marginal (e.g., upper) diameters. For cases where tree heights are not measured, methods for estimating additional variance in volume predictions resulting from uncertainties in tree height models (tariffs) are provided. The example data include the taper curve parameters for Norway spruce used in the 3rd German NFI fitted to 380 trees and a subset of section-wise diameter measurements of these trees. The functions implemented here are detailed in Kublin, E., Breidenbach, J., Kaendler, G. (2013) <doi:10.1007/s10342-013-0715-0>. |
Authors: | Edgar Kublin [aut], Johannes Breidenbach [ctb], Christian Vonderach [ctb, cre] |
Maintainer: | Christian Vonderach <[email protected]> |
License: | GPL (>= 2) |
Version: | 0.5.3 |
Built: | 2025-01-24 03:24:23 UTC |
Source: | https://github.com/vochr/taper |
This package implements functions for fitting taper curves (a semiparametric linear mixed effects taper model) to diameter measurements along stems. Further functions are provided to estimate the variance/confidence intervals around the predicted curves, to calculate timber volume (also by sections) and marginal (e.g., upper) diameters. For cases where tree heights are not measured, methods for estimating additional variance in volume predictions resulting from uncertainties in tree height models (tariffs) are provided. The example data include the taper curve parameters for Norway spruce used in the 3rd German NFI fitted to 380 trees and a subset of section-wise diameter measurements of these trees.
Package: | TapeR |
Type: | Package |
License: | GPL-2 |
LazyLoad: | yes |
Fits taper models using diameter measurements along the stem. Uses fitted models and arbitrary numbers of diameter measurements to estimate diameter at any position along the stem. Estimates timber volume from the taper curve. Provides variances for all estimates.
Edgar Kublin
Maintainer: Johannes Breidenbach <[email protected]>
Kublin, E., Breidenbach, J., Kaendler, G. (2013) A flexible stem taper and volume prediction method based on mixed-effects B-spline regression, Eur J For Res, 132:983-997.
TapeR_FIT_LME.f
, E_DHx_HmDm_HT.f
,
DxHx.df
, SK.par.lme
, HT.par
Internal function not usually called by users
BSplines(knots = c(seq(0, 1, 0.1)), ord = 4, der = 0, x = c(seq(0, 1, 0.01)), ...)
BSplines(knots = c(seq(0, 1, 0.1)), ord = 4, der = 0, x = c(seq(0, 1, 0.01)), ...)
knots |
knot positions for spline function |
ord |
order of the spline function |
der |
derivatives |
x |
height measurements |
... |
not currently used |
internally splineDesign
is called
B-Splines matrix build using splineDesign
Edgar Kublin
Internal function not usually called by users
CdN_DHxHt.f(Ht, Hx, qD, Hm, Dm, par.lme, Rfn, ...)
CdN_DHxHt.f(Ht, Hx, qD, Hm, Dm, par.lme, Rfn, ...)
Ht |
tree height |
Hx |
Numeric vector of stem heights (m) along which to return the expected diameter. |
qD |
vector of quantiles, passed to |
Hm |
measured height of respective diameters |
Dm |
measured diameter |
par.lme |
List of taper model parameters obtained by |
Rfn |
list with function name to provide estimated or assumed residual variances for the given measurements, optionally parameters for such functions |
... |
not currently used |
percentile for estimated taper curve diameter at position Hx
given Ht
, Hm
and Dm
Edgar Kublin
Internal function not usually called by users
dN.f(x, mw, sd, ...)
dN.f(x, mw, sd, ...)
x |
vector of quantiles |
mw |
vector of means |
sd |
vector of standard deviations |
... |
not currently used |
numeric density of normal distribution
Edgar Kublin
Example dataset of 10 trees with 10 diameter and height measurements for each tree.
data(DxHx.df)
data(DxHx.df)
A data frame with 172 observations on the following 4 variables.
Id
Numeric vector of tree IDs.
Dx
Numeric vector of diameter measurements.
Hx
Numeric vector of height measurements.
Ht
Numeric vector of tree height (repeated for each measurement in each tree).
Measured for BWI3.
Kublin, E., Breidenbach, J., Kaendler, G. (2013) A flexible stem taper and volume prediction method based on mixed-effects B-spline regression, Eur J For Res, 132:983-997.
data(DxHx.df) head(DxHx.df)
data(DxHx.df) head(DxHx.df)
Calibrates a taper curve based on at least one diameter measurement and returns the expected diameters and exact variances
E_DHx_HmDm_HT_CIdHt.f( Hx, Hm, Dm, mHt, sHt, par.lme, Rfn = list(fn = "sig2"), ... )
E_DHx_HmDm_HT_CIdHt.f( Hx, Hm, Dm, mHt, sHt, par.lme, Rfn = list(fn = "sig2"), ... )
Hx |
Numeric vector of stem heights (m) along which to return the expected diameter. |
Hm |
Numeric vector of stem heights (m) along which diameter
measurements were taken for calibration. Can be of length 1. Must be of same
length as |
Dm |
Numeric vector of diameter measurements (cm) taken for calibration.
Can be of length 1. Must be of same length as |
mHt |
Scalar. Tree height (m). |
sHt |
Scalar. Standard deviation of stem height. Can be 0 if height was measured without error. |
par.lme |
List of taper model parameters obtained by
|
Rfn |
list with function name to provide estimated or assumed residual variances for the given measurements, optionally parameters for such functions |
... |
not currently used |
calibrates the tree specific taper curve and calculates 'exact' confidence intervals, which can be useful for plotting. Attention: this function is somewhat time-consuming.
a matrix with six columns:
Hx: Numeric vector of heights (m) along which to return the expected diameter.
q_DHx_u: Lower confidence interval (cm). (95% CI except for estimates close to the stem tip.)
DHx: Diameter estimate (cm).
q_DHx_o: Upper CI (cm).
cP_DHx_u: Probability of observations <q_DHx_u
.
cP_DHx_o: Probability of observations <q_DHx_o
.
Edgar Kublin
Kublin, E., Breidenbach, J., Kaendler, G. (2013) A flexible stem taper and volume prediction method based on mixed-effects B-spline regression, Eur J For Res, 132:983-997.
# example data data(DxHx.df) # taper curve parameters based on all measured trees data(SK.par.lme) #select data of first tree Idi <- (DxHx.df[,"Id"] == unique(DxHx.df$Id)[1]) (tree1 <- DxHx.df[Idi,]) ## Predict the taper curve based on the diameter measurement in 2 m ## height and known height tc.tree1 <- E_DHx_HmDm_HT.f(Hx=1:tree1$Ht[1], Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 0, par.lme = SK.par.lme) #plot the predicted taper curve plot(tc.tree1$Hx, tc.tree1$DHx, type="l", las=1) #lower CI lines(tc.tree1$Hx, tc.tree1$CI_Mean[,1], lty=2) #upper CI lines(tc.tree1$Hx, tc.tree1$CI_Mean[,3], lty=2) #lower prediction interval lines(tc.tree1$Hx, tc.tree1$CI_Pred[,1], lty=3) #upper prediction interval lines(tc.tree1$Hx, tc.tree1$CI_Pred[,3], lty=3) #add measured diameter used for calibration points(tree1$Hx[3], tree1$Dx[3], pch=3, col=2) #add the observations points(tree1$Hx, tree1$Dx) ## Calculate "exact" CIs. Careful: This takes a while! #library(pracma)# for numerical integration with gaussLegendre() tc.tree1.exact <- E_DHx_HmDm_HT_CIdHt.f(Hx=1:tree1$Ht[1], Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt=tree1$Ht[1], sHt=1, par.lme=SK.par.lme) #add exact confidence intervals to approximate intervals above - fits #quite well lines(tc.tree1.exact[,1], tc.tree1.exact[,2], lty=2,col=2) lines(tc.tree1.exact[,1], tc.tree1.exact[,4], lty=2,col=2)
# example data data(DxHx.df) # taper curve parameters based on all measured trees data(SK.par.lme) #select data of first tree Idi <- (DxHx.df[,"Id"] == unique(DxHx.df$Id)[1]) (tree1 <- DxHx.df[Idi,]) ## Predict the taper curve based on the diameter measurement in 2 m ## height and known height tc.tree1 <- E_DHx_HmDm_HT.f(Hx=1:tree1$Ht[1], Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 0, par.lme = SK.par.lme) #plot the predicted taper curve plot(tc.tree1$Hx, tc.tree1$DHx, type="l", las=1) #lower CI lines(tc.tree1$Hx, tc.tree1$CI_Mean[,1], lty=2) #upper CI lines(tc.tree1$Hx, tc.tree1$CI_Mean[,3], lty=2) #lower prediction interval lines(tc.tree1$Hx, tc.tree1$CI_Pred[,1], lty=3) #upper prediction interval lines(tc.tree1$Hx, tc.tree1$CI_Pred[,3], lty=3) #add measured diameter used for calibration points(tree1$Hx[3], tree1$Dx[3], pch=3, col=2) #add the observations points(tree1$Hx, tree1$Dx) ## Calculate "exact" CIs. Careful: This takes a while! #library(pracma)# for numerical integration with gaussLegendre() tc.tree1.exact <- E_DHx_HmDm_HT_CIdHt.f(Hx=1:tree1$Ht[1], Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt=tree1$Ht[1], sHt=1, par.lme=SK.par.lme) #add exact confidence intervals to approximate intervals above - fits #quite well lines(tc.tree1.exact[,1], tc.tree1.exact[,2], lty=2,col=2) lines(tc.tree1.exact[,1], tc.tree1.exact[,4], lty=2,col=2)
Calibrates a taper curve based on at least one diameter measurement and returns the expected diameters and approximate variances
E_DHx_HmDm_HT.f( Hx, Hm, Dm, mHt, sHt = 0, par.lme, Rfn = list(fn = "sig2"), ... )
E_DHx_HmDm_HT.f( Hx, Hm, Dm, mHt, sHt = 0, par.lme, Rfn = list(fn = "sig2"), ... )
Hx |
Numeric vector of stem heights (m) along which to return the expected diameter. |
Hm |
Numeric vector of stem heights (m) along which diameter
measurements were taken for calibration. Can be of length 1. Must be of same
length as |
Dm |
Numeric vector of diameter measurements (cm) taken for calibration.
Can be of length 1. Must be of same length as |
mHt |
Scalar. Tree height (m). |
sHt |
Scalar. Standard deviation of stem height. Can be 0 if height was measured without error. |
par.lme |
List of taper model parameters obtained by
|
Rfn |
list with function name to provide estimated or assumed residual variances for the given measurements, optionally parameters for such functions |
... |
not currently used |
calibrates the tree specific taper curve and calculates approximate
confidence intervals, which can be useful for plotting. Uncertainty resulting
from tariff height estimates if tree height was not measured is incorporated.
Using Rfn
the taper curve can be forced through the measured
diameters, c.f. resVar
.
a list holding nine elements:
DHx: Numeric vector of diameters (cm) (expected value) along the
heights given by Hx
.
Hx: Numeric vector of heights (m) along which to return the expected diameter.
MSE_Mean: Mean squared error for the expected value of the diameter.
CI_Mean: Confidence interval. Matrix of the 95% conf. int. for the expected value of the diameter (cm). First column: lower limit, second column: mean, third column: upper limit.
KOV_Mean: Variance-Covariance matrix for the expected value of the diameter.
MSE_Pred: Mean squared error for the prediction of the diameter.
CI_Pred: Prediction interval. Matrix of the 95% conf. int. for the prediction of the diameter (cm). First column: lower limit, second column: mean, third column: upper limit.
KOV_Pred: Variance-Covariance matrix for the prediction of the diameter.
Rfn: Function applied for estimated or assumed residual variance.
Edgar Kublin and Christian Vonderach
Kublin, E., Breidenbach, J., Kaendler, G. (2013) A flexible stem taper and volume prediction method based on mixed-effects B-spline regression, Eur J For Res, 132:983-997.
# example data data(DxHx.df) # taper curve parameters based on all measured trees data(SK.par.lme) #select data of first tree Idi <- (DxHx.df[,"Id"] == unique(DxHx.df$Id)[1]) (tree1 <- DxHx.df[Idi,]) ## Predict the taper curve based on the diameter measurement in 2 m ## height and known height tc.tree1 <- E_DHx_HmDm_HT.f(Hx=1:tree1$Ht[1], Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 0, par.lme = SK.par.lme) #plot the predicted taper curve plot(tc.tree1$Hx, tc.tree1$DHx, type="l", las=1) #lower CI lines(tc.tree1$Hx, tc.tree1$CI_Mean[,1], lty=2) #upper CI lines(tc.tree1$Hx, tc.tree1$CI_Mean[,3], lty=2) #lower prediction interval lines(tc.tree1$Hx, tc.tree1$CI_Pred[,1], lty=3) #upper prediction interval lines(tc.tree1$Hx, tc.tree1$CI_Pred[,3], lty=3) #add measured diameter used for calibration points(tree1$Hx[3], tree1$Dx[3], pch=3, col=2) #add the observations points(tree1$Hx, tree1$Dx) ## feature of forcing taper curve through measured diameters i <- c(3, 6) tc.tree1 <- E_DHx_HmDm_HT.f(Hx=seq(0, tree1$Ht[1], 0.1), Hm=tree1$Hx[i], Dm=tree1$Dx[i], mHt = tree1$Ht[1], sHt = 0, par.lme = SK.par.lme, Rfn=list(fn="sig2")) tc.tree2 <- E_DHx_HmDm_HT.f(Hx=seq(0, tree1$Ht[1], 0.1), Hm=tree1$Hx[i], Dm=tree1$Dx[i], mHt = tree1$Ht[1], sHt = 0, par.lme = SK.par.lme, Rfn=list(fn="zero")) # plot the predicted taper curve plot(tc.tree1$Hx, tc.tree1$DHx, type="l", las=1) # added taper curve through measurement points(x=tc.tree2$Hx, y=tc.tree2$DHx, type="l", lty=2) # closer window plot(tc.tree1$Hx, tc.tree1$DHx, type="l", las=1, xlim=c(0, 8), ylim=c(24, 30)) # added taper curve through measurement points(x=tc.tree2$Hx, y=tc.tree2$DHx, type="l", lty=2) # add measured diameter used for calibration points(tree1$Hx[i], tree1$Dx[i], pch=3, col=2) # add the observations points(tree1$Hx, tree1$Dx) ## apply yet another residual variance function i <- c(1, 2, 3) # calibrating with 0.5, 1m and 2m, assuming no error in 0.5m zrv <- tree1$Hx[1] / tree1$Ht[1] # assumed zero resiudal variance # assumed residual variance per measurement TapeR:::resVar(relH = tree1$Hx[i] / tree1$Ht[1], fn = "dlnorm", sig2 = SK.par.lme$sig2_eps, par = list(zrv=zrv)) tc.tree3 <- E_DHx_HmDm_HT.f(Hx=seq(0, tree1$Ht[1], 0.1), Hm=tree1$Hx[i], Dm=tree1$Dx[i], mHt = tree1$Ht[1], sHt = 0, par.lme = SK.par.lme, Rfn=list(fn="dlnorm", par=list(zrv=zrv))) plot(tc.tree1$Hx, tc.tree1$DHx, type="l", las=1, xlim=c(0, 4)) points(x=tc.tree3$Hx, y=tc.tree3$DHx, type="l", lty=2) points(tree1$Hx[i], tree1$Dx[i], pch=3, col=2) points(tree1$Hx, tree1$Dx)
# example data data(DxHx.df) # taper curve parameters based on all measured trees data(SK.par.lme) #select data of first tree Idi <- (DxHx.df[,"Id"] == unique(DxHx.df$Id)[1]) (tree1 <- DxHx.df[Idi,]) ## Predict the taper curve based on the diameter measurement in 2 m ## height and known height tc.tree1 <- E_DHx_HmDm_HT.f(Hx=1:tree1$Ht[1], Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 0, par.lme = SK.par.lme) #plot the predicted taper curve plot(tc.tree1$Hx, tc.tree1$DHx, type="l", las=1) #lower CI lines(tc.tree1$Hx, tc.tree1$CI_Mean[,1], lty=2) #upper CI lines(tc.tree1$Hx, tc.tree1$CI_Mean[,3], lty=2) #lower prediction interval lines(tc.tree1$Hx, tc.tree1$CI_Pred[,1], lty=3) #upper prediction interval lines(tc.tree1$Hx, tc.tree1$CI_Pred[,3], lty=3) #add measured diameter used for calibration points(tree1$Hx[3], tree1$Dx[3], pch=3, col=2) #add the observations points(tree1$Hx, tree1$Dx) ## feature of forcing taper curve through measured diameters i <- c(3, 6) tc.tree1 <- E_DHx_HmDm_HT.f(Hx=seq(0, tree1$Ht[1], 0.1), Hm=tree1$Hx[i], Dm=tree1$Dx[i], mHt = tree1$Ht[1], sHt = 0, par.lme = SK.par.lme, Rfn=list(fn="sig2")) tc.tree2 <- E_DHx_HmDm_HT.f(Hx=seq(0, tree1$Ht[1], 0.1), Hm=tree1$Hx[i], Dm=tree1$Dx[i], mHt = tree1$Ht[1], sHt = 0, par.lme = SK.par.lme, Rfn=list(fn="zero")) # plot the predicted taper curve plot(tc.tree1$Hx, tc.tree1$DHx, type="l", las=1) # added taper curve through measurement points(x=tc.tree2$Hx, y=tc.tree2$DHx, type="l", lty=2) # closer window plot(tc.tree1$Hx, tc.tree1$DHx, type="l", las=1, xlim=c(0, 8), ylim=c(24, 30)) # added taper curve through measurement points(x=tc.tree2$Hx, y=tc.tree2$DHx, type="l", lty=2) # add measured diameter used for calibration points(tree1$Hx[i], tree1$Dx[i], pch=3, col=2) # add the observations points(tree1$Hx, tree1$Dx) ## apply yet another residual variance function i <- c(1, 2, 3) # calibrating with 0.5, 1m and 2m, assuming no error in 0.5m zrv <- tree1$Hx[1] / tree1$Ht[1] # assumed zero resiudal variance # assumed residual variance per measurement TapeR:::resVar(relH = tree1$Hx[i] / tree1$Ht[1], fn = "dlnorm", sig2 = SK.par.lme$sig2_eps, par = list(zrv=zrv)) tc.tree3 <- E_DHx_HmDm_HT.f(Hx=seq(0, tree1$Ht[1], 0.1), Hm=tree1$Hx[i], Dm=tree1$Dx[i], mHt = tree1$Ht[1], sHt = 0, par.lme = SK.par.lme, Rfn=list(fn="dlnorm", par=list(zrv=zrv))) plot(tc.tree1$Hx, tc.tree1$DHx, type="l", las=1, xlim=c(0, 4)) points(x=tc.tree3$Hx, y=tc.tree3$DHx, type="l", lty=2) points(tree1$Hx[i], tree1$Dx[i], pch=3, col=2) points(tree1$Hx, tree1$Dx)
Calibrates a taper curve based on at least one diameter measurement and returns the height of a given diameter
E_HDx_HmDm_HT.f( Dx, Hm, Dm, mHt, sHt = 0, par.lme, Rfn = list(fn = "sig2"), ... )
E_HDx_HmDm_HT.f( Dx, Hm, Dm, mHt, sHt = 0, par.lme, Rfn = list(fn = "sig2"), ... )
Dx |
Scalar. Diameter for which to return height. |
Hm |
Numeric vector of stem heights (m) along which diameter
measurements were taken for calibration. Can be of length 1. Must be of same
length as |
Dm |
Numeric vector of diameter measurements (cm) taken for calibration.
Can be of length 1. Must be of same length as |
mHt |
Scalar. Tree height (m). |
sHt |
Scalar. Standard deviation of stem height. Can be 0 if height was measured without error. |
par.lme |
List of taper model parameters obtained by
|
Rfn |
list with function name to provide estimated or assumed residual variances for the given measurements, optionally parameters for such functions |
... |
not currently used |
returns the height given a certain diameter.
A scalar. Estimated height (m) given a diameter.
Edgar Kublin
Kublin, E., Breidenbach, J., Kaendler, G. (2013) A flexible stem taper and volume prediction method based on mixed-effects B-spline regression, Eur J For Res, 132:983-997.
# example data data(DxHx.df) # taper curve parameters based on all measured trees data(SK.par.lme) #select data of first tree Idi <- (DxHx.df[,"Id"] == unique(DxHx.df$Id)[1]) (tree1 <- DxHx.df[Idi,]) ## Predict the taper curve based on the diameter measurement in 2 m ## height and known height tc.tree1 <- E_DHx_HmDm_HT.f(Hx=1:tree1$Ht[1], Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 0, par.lme = SK.par.lme) #plot the predicted taper curve plot(tc.tree1$Hx, tc.tree1$DHx, type="l", las=1) ## Calculate the height given a certain diameter threshold, say 8.5 cm ht.tree1.d8.5 <- E_HDx_HmDm_HT.f (Dx=8.5, Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 1, par.lme = SK.par.lme, Rfn=list(fn="sig2")) # add to plot points(x=ht.tree1.d8.5, y=8.5, pch=8, col=2, cex=2)
# example data data(DxHx.df) # taper curve parameters based on all measured trees data(SK.par.lme) #select data of first tree Idi <- (DxHx.df[,"Id"] == unique(DxHx.df$Id)[1]) (tree1 <- DxHx.df[Idi,]) ## Predict the taper curve based on the diameter measurement in 2 m ## height and known height tc.tree1 <- E_DHx_HmDm_HT.f(Hx=1:tree1$Ht[1], Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 0, par.lme = SK.par.lme) #plot the predicted taper curve plot(tc.tree1$Hx, tc.tree1$DHx, type="l", las=1) ## Calculate the height given a certain diameter threshold, say 8.5 cm ht.tree1.d8.5 <- E_HDx_HmDm_HT.f (Dx=8.5, Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 1, par.lme = SK.par.lme, Rfn=list(fn="sig2")) # add to plot points(x=ht.tree1.d8.5, y=8.5, pch=8, col=2, cex=2)
Estimate volume for a complete stem from bottom to tip or for a section defined by lower and upper diameter or height. Variances for estimated volumes are calculated.
E_VOL_AB_HmDm_HT.f( Hm, Dm, mHt, sHt = 0, A = NULL, B = NULL, iDH = "D", par.lme, Rfn = list(fn = "sig2"), IA = F, nGL = 51, ... )
E_VOL_AB_HmDm_HT.f( Hm, Dm, mHt, sHt = 0, A = NULL, B = NULL, iDH = "D", par.lme, Rfn = list(fn = "sig2"), IA = F, nGL = 51, ... )
Hm |
Numeric vector of stem heights (m) along which diameter
measurements were taken for calibration. Can be of length 1. Must be of same
length as |
Dm |
Numeric vector of diameter measurements (cm) taken for calibration.
Can be of length 1. Must be of same length as |
mHt |
Scalar. Tree height (m). |
sHt |
Scalar. Standard deviation of stem height. Can be 0 if height was measured without error. |
A |
Numeric scalar defining the lower threshold of a stem section for
volume estimation. Depends on |
B |
Numeric scalar defining the upper threshold of a stem section for
volume estimation. Depends on |
iDH |
Character scalar. Either "D" or "H". Type of threshold for section
volume estimation. See |
par.lme |
List of taper model parameters obtained by
|
Rfn |
list with function name to provide estimated or assumed residual variances for the given measurements, optionally parameters for such functions |
IA |
Logic scalar. If TRUE, variance calculation of height estimate based on 2-point distribution. If FALSE, variance calculation of height estimate based on Normal approximation. |
nGL |
Numeric scalar. Number of support points for numerical integration. |
... |
not currently used |
calculates the volume for a complete stem or sections defined by
A
and B
, which might be defined as diameter or height. The
parameter Rfn
can be used to force the taper curve through the
measured points (e.g. by Rfn=list(fn="zero")
, cf. resVar
).
a list holding nine elements:
E_VOL: Estimated volume (m^3).
VAR_VOL: Variance of the volume estimate.
Hm: Height of diameter measurement (m).
Dm: Diameter measurement (cm).
Ht: Tree height (m).
Da: Diameter at lower section threshold (cm).
Db: Diameter at upper section threshold (cm).
Ha: Height at lower section threshold (m).
Hb: Height at upper section threshold (m).
Rfn: Function applied for estimated or assumed residual variance.
Edgar Kublin
Kublin, E., Breidenbach, J., Kaendler, G. (2013) A flexible stem taper and volume prediction method based on mixed-effects B-spline regression, Eur J For Res, 132:983-997.
# example data data(DxHx.df) # taper curve parameters based on all measured trees data(SK.par.lme) #select data of first tree Idi <- (DxHx.df[,"Id"] == unique(DxHx.df$Id)[1]) (tree1 <- DxHx.df[Idi,]) ## Calculate the timber volume for the whole stem VOL <- E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 0, # no height variance assumed par.lme = SK.par.lme) VOL$E_VOL #' expected value VOL$VAR_VOL #' corresponding variance VOL$Rfn ## Calculate the timber volume for the whole stem, using Rfn="zero" VOL <- E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 0, # no height variance assumed par.lme = SK.par.lme, Rfn = list(fn="zero")) VOL$E_VOL #' expected value VOL$VAR_VOL #' corresponding variance VOL$Rfn ## Calculate the timber volume for the whole stem VOL <- E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 1, # no height variance assumed par.lme = SK.par.lme) VOL$E_VOL #' expected value VOL$VAR_VOL #' corresponding variance VOL$Rfn ## Calculate the timber volume for the whole stem, using Rfn="zero" VOL <- E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 1, # height variance assumed par.lme = SK.par.lme, Rfn = list(fn="zero")) VOL$E_VOL #' expected value VOL$VAR_VOL #' corresponding variance VOL$Rfn ## Calculate the timber volume for a selected section given a height (0.3 - 5 m) VOL <- E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 1, par.lme = SK.par.lme, A=0.3, B=5, iDH = "H") VOL$E_VOL #' expected value VOL$VAR_VOL #' corresponding variance VOL$Rfn ## Calculate the timber volume for a selected section given a height (0.3 - 5 m) VOL <- E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 1, par.lme = SK.par.lme, A=0.3, B=5, iDH = "H", Rfn=list(fn="zero")) VOL$E_VOL #' expected value VOL$VAR_VOL #' corresponding variance VOL$Rfn ## Calculate the timber volume for a selected section given a diameter ## threshold (30cm - 15cm) (negative value if A<B) VOL <- E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 1, par.lme = SK.par.lme, A=30, B=15, iDH = "D") VOL$E_VOL #' expected value VOL$VAR_VOL #' corresponding variance ## Not run: ## The variance estimate resulting from the tree height uncertainty using ## a Normal approximation takes much longer... ptm <- proc.time() E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 1, par.lme = SK.par.lme, IA=FALSE) proc.time() - ptm ##... than the calculation using a 2-point distribution... ptm <- proc.time() E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 1, par.lme = SK.par.lme, IA=TRUE) proc.time() - ptm ##...fastest if no height variance is assumed ptm <- proc.time() E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 0, par.lme = SK.par.lme, IA=FALSE) proc.time() - ptm ## Also the number of supportive points for the numerical integration ## influences the calculation time ptm <- proc.time() E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 0, par.lme = SK.par.lme, IA=FALSE, nGL=10) proc.time() - ptm ##' End(Not run)
# example data data(DxHx.df) # taper curve parameters based on all measured trees data(SK.par.lme) #select data of first tree Idi <- (DxHx.df[,"Id"] == unique(DxHx.df$Id)[1]) (tree1 <- DxHx.df[Idi,]) ## Calculate the timber volume for the whole stem VOL <- E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 0, # no height variance assumed par.lme = SK.par.lme) VOL$E_VOL #' expected value VOL$VAR_VOL #' corresponding variance VOL$Rfn ## Calculate the timber volume for the whole stem, using Rfn="zero" VOL <- E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 0, # no height variance assumed par.lme = SK.par.lme, Rfn = list(fn="zero")) VOL$E_VOL #' expected value VOL$VAR_VOL #' corresponding variance VOL$Rfn ## Calculate the timber volume for the whole stem VOL <- E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 1, # no height variance assumed par.lme = SK.par.lme) VOL$E_VOL #' expected value VOL$VAR_VOL #' corresponding variance VOL$Rfn ## Calculate the timber volume for the whole stem, using Rfn="zero" VOL <- E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 1, # height variance assumed par.lme = SK.par.lme, Rfn = list(fn="zero")) VOL$E_VOL #' expected value VOL$VAR_VOL #' corresponding variance VOL$Rfn ## Calculate the timber volume for a selected section given a height (0.3 - 5 m) VOL <- E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 1, par.lme = SK.par.lme, A=0.3, B=5, iDH = "H") VOL$E_VOL #' expected value VOL$VAR_VOL #' corresponding variance VOL$Rfn ## Calculate the timber volume for a selected section given a height (0.3 - 5 m) VOL <- E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 1, par.lme = SK.par.lme, A=0.3, B=5, iDH = "H", Rfn=list(fn="zero")) VOL$E_VOL #' expected value VOL$VAR_VOL #' corresponding variance VOL$Rfn ## Calculate the timber volume for a selected section given a diameter ## threshold (30cm - 15cm) (negative value if A<B) VOL <- E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 1, par.lme = SK.par.lme, A=30, B=15, iDH = "D") VOL$E_VOL #' expected value VOL$VAR_VOL #' corresponding variance ## Not run: ## The variance estimate resulting from the tree height uncertainty using ## a Normal approximation takes much longer... ptm <- proc.time() E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 1, par.lme = SK.par.lme, IA=FALSE) proc.time() - ptm ##... than the calculation using a 2-point distribution... ptm <- proc.time() E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 1, par.lme = SK.par.lme, IA=TRUE) proc.time() - ptm ##...fastest if no height variance is assumed ptm <- proc.time() E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 0, par.lme = SK.par.lme, IA=FALSE) proc.time() - ptm ## Also the number of supportive points for the numerical integration ## influences the calculation time ptm <- proc.time() E_VOL_AB_HmDm_HT.f(Hm=tree1$Hx[3], Dm=tree1$Dx[3], mHt = tree1$Ht[1], sHt = 0, par.lme = SK.par.lme, IA=FALSE, nGL=10) proc.time() - ptm ##' End(Not run)
Internal function not usually called by users
EYx_ssp.f(knt, coe, x, ...)
EYx_ssp.f(knt, coe, x, ...)
knt |
knots position of B-Splines |
coe |
estimated coefficient for B-Splines |
x |
position at which to evaluate B-Splines model |
... |
not currently used |
expected diameter given knots and coefficients at position (height)
x
.
Edgar Kublin
Height is only measured on a subset of the trees on a sample plots. This Height tariff is used to estimate the height of the trees with only a dbh measurement.
data(HT.par)
data(HT.par)
The format is: List of 4 $ knt.mw: num [1:16] 0 0 0 0 19.6 ... $ coe.mw: num [1:12] 1.3 7.28 15.1 21.75 24.39 ... $ knt.sd: num [1:53] 0 0 0 0 7.52 ... $ coe.sd: num [1:49] 0 0.618 1.376 2.142 2.486 ...
Kublin, E., Breidenbach, J., Kaendler, G. (2013) A flexible stem taper and volume prediction method based on mixed-effects B-spline regression, Eur J For Res, 132:983-997.
data(HT.par) ## maybe str(HT.par) ; plot(HT.par) ...
data(HT.par) ## maybe str(HT.par) ; plot(HT.par) ...
Internal function not usually called by users
Hx_root.f(Hx, Dx, Hm, Dm, mHt, sHt, par.lme, Rfn, ...)
Hx_root.f(Hx, Dx, Hm, Dm, mHt, sHt, par.lme, Rfn, ...)
Hx |
Numeric vector of stem heights (m) along which to return the expected diameter |
Dx |
expected diameter |
Hm |
Numeric vector of stem heights (m) along which diameter measurements
were taken for calibration. Can be of length 1. Must be of same length as |
Dm |
Numeric vector of diameter measurements (cm) taken for calibration.
Can be of length 1. Must be of same length as |
mHt |
Scalar. Tree height (m) |
sHt |
Scalar. Standard deviation of stem height. Can be 0 if height was measured without error |
par.lme |
List of taper model parameters obtained by
|
Rfn |
list with function name to provide estimated or assumed residual variances for the given measurements, optionally parameters for such functions |
... |
not currently used |
function is called by uniroot
inside
E_HDx_HmDm_HT.f
deviation between observed diameter Dx
and diameter in height
Hx
.
Edgar Kublin
Internal function not usually called by users
Int_CdN_DHx_dHt.f(qD, Hx, Hm, Dm, mHt, sHt, par.lme, Rfn, nGL = 51, ...)
Int_CdN_DHx_dHt.f(qD, Hx, Hm, Dm, mHt, sHt, par.lme, Rfn, nGL = 51, ...)
qD |
vector of quantiles, finally passed to |
Hx |
Numeric vector of stem heights (m) along which to return the expected diameter |
Hm |
Numeric vector of stem heights (m) along which diameter measurements
were taken for calibration. Can be of length 1. Must be of same length as |
Dm |
Numeric vector of diameter measurements (cm) taken for calibration.
Can be of length 1. Must be of same length as |
mHt |
Scalar. Tree height (m) |
sHt |
Scalar. Standard deviation of stem height. Can be 0 if height was measured without error |
par.lme |
List of taper model parameters obtained by |
Rfn |
list with function name to provide estimated or assumed residual variances for the given measurements, optionally parameters for such functions |
nGL |
Numeric scalar. Number of support points for numerical integration |
... |
not currently used |
Int_CdN_dN
Edgar Kublin
Internal function not usually called by users
Int_E_VOL_AB_HmDm_HT_dHt.f( Hm, Dm, A = NULL, B = NULL, iDH = "D", mw_HtT, sd_HtT, par.lme, Rfn = list(fn = "sig2"), IA = FALSE, nGL = 51, ... )
Int_E_VOL_AB_HmDm_HT_dHt.f( Hm, Dm, A = NULL, B = NULL, iDH = "D", mw_HtT, sd_HtT, par.lme, Rfn = list(fn = "sig2"), IA = FALSE, nGL = 51, ... )
Hm |
Numeric vector of stem heights (m) along which diameter measurements
were taken for calibration. Can be of length 1. Must be of same length as |
Dm |
Numeric vector of diameter measurements (cm) taken for calibration.
Can be of length 1. Must be of same length as |
A |
Numeric scalar defining the lower threshold of a stem section for volume
estimation. Depends on |
B |
Numeric scalar defining the upper threshold of a stem section for volume
estimation. Depends on |
iDH |
Character scalar. Either "D" or "H". Type of threshold for
section volume estimation. See |
mw_HtT |
Scalar. Tree height (m) |
sd_HtT |
Scalar. Standard deviation of stem height. Can be 0 if height was measured without error |
par.lme |
List of taper model parameters obtained by |
Rfn |
list with function name to provide estimated or assumed residual variances for the given measurements, optionally parameters for such functions |
IA |
Logic scalar. If TRUE, variance calculation of height estimate based on 2-point distribution. If FALSE, variance calculation of height estimate based on Normal approximation. |
nGL |
Numeric scalar. Number of support points for numerical integration. |
... |
not currently used |
integrating the taper curve considering uncertainty of height measurement
list with expected volume, variance of volume and squared expected value incorporating the uncertainty of height measurement
Edgar Kublin
Internal function not usually called by users
qD.rout.f(qD, alpha = 0.975, Hx, Hm, Dm, mHt, sHt, par.lme, Rfn, nGL = 51, ...)
qD.rout.f(qD, alpha = 0.975, Hx, Hm, Dm, mHt, sHt, par.lme, Rfn, nGL = 51, ...)
qD |
vector of quantiles, finally passed to |
alpha |
quantile for which root is sought |
Hx |
Numeric vector of stem heights (m) along which to return the expected diameter |
Hm |
Numeric vector of stem heights (m) along which diameter measurements
were taken for calibration. Can be of length 1. Must be of same length as |
Dm |
Numeric vector of diameter measurements (cm) taken for calibration.
Can be of length 1. Must be of same length as |
mHt |
Scalar. Tree height (m) |
sHt |
Scalar. Standard deviation of stem height. Can be 0 if height was measured without error |
par.lme |
List of taper model parameters obtained by
|
Rfn |
list with function name to provide estimated or assumed residual variances for the given measurements, optionally parameters for such functions |
nGL |
Numeric scalar. Number of support points for numerical integration. |
... |
not currently used |
qD
for given alpha
with respect to Int_CdN_DHx_dHt.f
Edgar Kublin
When estimating a tree specific taper curve based on given measurements one can modify the assumed measurement uncertainty by these functions
resVar(relH, fn, sig2, par = NULL)
resVar(relH, fn, sig2, par = NULL)
relH |
relative tree height for which the assumed residual variance should be calculated |
fn |
name of function to be applied as character string |
sig2 |
residual variance from fitted model, cf.
|
par |
either NULL or a list with parameters to the different functions. See details. |
When estimating a tree specific taper curve based on given measurements the residual variance of the model is taken into account to estimate the tree specific random effects. Alternatively, it is possible to make assumptions about the measurement error, eventually at specific relative heights. With that, one can e.g. force the taper curve through the given measurements. Standard behaviour not necessarily leads to passing the measurements, if more than one is given.
Different functions are available. sig2
applies the model residual
variance and hence is the default behaviour.
zero
means assuming no residual variance and forcing the taper curve
through the given measurements. Care has to be taken in this case because
forcing the taper curve through a lot of measurements might result in
implausible results.
linear
interpolates between zero and the given residual variance along
the stem, i.e. from bottom to tree top.
bilinear
puts zero variance not at zero but at a predefined location
(can be given via par
). Below and above a linear interpolation is
done up to the given residual variance. If zero variance position is not
given, it is set at 5% of tree height (approximately height of dbh).
laglinear
assumes zero variance up to a predefined location (defaults
to 5% of tree height) and interpolates upwards to the given residual
variance of the model.
quadratic
function distributes residual variance according to a
quadratic function along the stem. It is build so that zero variance is put
at a predefined location (defaults to 5% of tree height) and model residual
variance (as a default) at tree top.
dnorm
and dlnorm
put residual variance in form of an inverse
normal or an inverse log-normal distribution along the stem with a
zero-minimum at a predefined location (defaults to 5% of tree height).
See examples for a visualisation.
For all functions (except zero
and sig2
) the point of zero
residual variance is defined by par$zrv
if given, otherwise set to 0.05.
For dnorm
one can additionally provide the parameter sd
to
determine standard deviation. By default it is set to zrv/3
; in case
of dlnorm
one can define lsd
(sdlog, cf. dlnorm
),
which is by default set to 1-sqrt(zrv)
. It is up to the user to define
meaningful parameters and use the functions in appropriate context.
vector of assumed residual variance
Christian Vonderach
curve(resVar(relH=x, fn = "sig2", 0.5)) curve(resVar(relH=x, fn = "zero", 0.5)) curve(resVar(relH=x, fn = "linear", 0.5)) curve(resVar(relH=x, fn = "bilinear", 0.5)) curve(resVar(relH=x, fn = "laglinear", 0.5)) curve(resVar(relH=x, fn = "quadratic", 0.5)) curve(resVar(relH=x, fn = "dnorm", 0.5)) curve(resVar(relH=x, fn = "dnorm", 0.5, par=list(zrv=0.2, sd=0.2/3))) curve(resVar(relH=x, fn = "dlnorm", 0.5)) curve(resVar(relH=x, fn = "dlnorm", 0.5, par=list(zrv=0.2))) invisible(sapply(seq(0.01, 0.99, length.out=20), function(a){ curve(resVar(relH=x, fn = "dlnorm", 0.5, par=list(zrv=a, lsd=(1-sqrt(a)))), n=1000) }))
curve(resVar(relH=x, fn = "sig2", 0.5)) curve(resVar(relH=x, fn = "zero", 0.5)) curve(resVar(relH=x, fn = "linear", 0.5)) curve(resVar(relH=x, fn = "bilinear", 0.5)) curve(resVar(relH=x, fn = "laglinear", 0.5)) curve(resVar(relH=x, fn = "quadratic", 0.5)) curve(resVar(relH=x, fn = "dnorm", 0.5)) curve(resVar(relH=x, fn = "dnorm", 0.5, par=list(zrv=0.2, sd=0.2/3))) curve(resVar(relH=x, fn = "dlnorm", 0.5)) curve(resVar(relH=x, fn = "dlnorm", 0.5, par=list(zrv=0.2))) invisible(sapply(seq(0.01, 0.99, length.out=20), function(a){ curve(resVar(relH=x, fn = "dlnorm", 0.5, par=list(zrv=a, lsd=(1-sqrt(a)))), n=1000) }))
This is the actual function to estimate diameters according to the fitted mixed B-splines model.
SK_EBLUP_LME.f(xm, ym, xp, par.lme, Rfn = list(fn = "sig2"), ...)
SK_EBLUP_LME.f(xm, ym, xp, par.lme, Rfn = list(fn = "sig2"), ...)
xm |
relative heights for which measurements are available |
ym |
corresponding diameter measurements in height |
xp |
relative heights for which predictions are required |
par.lme |
Fitted model object, return of |
Rfn |
list with function name to provide estimated or assumed residual variances for the given measurements, optionally parameters for such functions |
... |
not currently used |
This function is the actual working horse for prediction using the
fitted taper model. Based on the model par.lme
and the measured
diameters ym
and corresponding (relative) heights xm
of a
specific tree (there might be just one measurement), the random
effect parameters and subsequently diameters are estimated. Depending on the
parameter Rfn
, the calibrated taper curve is forced through the
given diameter ym
(Rfn = list(fn="zero")
), or calibrated using
the complete residual variance-covariance information
(Rfn = list(fn="sig2")
, the default).
Further assumptions are possible, see also resVar
and
Kublin et al. (2013) p. 987 for more details.
a list holding nine elements:
b_fix fixed effects parameter of taper model
b_rnd random effects parameter given tree (posterior mean b_k)
yp estimated diameter in height xp
KOV_Mean variance-covariance matrix of expected value
KOV_Pred variance-covariance matrix of prediction
CI_Mean mean and limits of confidence interval
MSE_Mean mean squared error of expected value
MSE_Pred mean squared error of prediction
CI_Pred mean and limits of prediction interval
Edgar Kublin
E_DHx_HmDm_HT.f
, E_VOL_AB_HmDm_HT.f
,
resVar
data("SK.par.lme") TapeR:::SK_EBLUP_LME.f(1.3/27, 30, 1.3/27, SK.par.lme) ## using empirical best linear unbiased estimator: estimate != 30 TapeR:::SK_EBLUP_LME.f(1.3/27, 30, 1.3/27, SK.par.lme, Rfn=list(fn="sig2"))$yp ## interpolate / force through given diameter: estimate == 30 TapeR:::SK_EBLUP_LME.f(1.3/27, 30, 1.3/27, SK.par.lme, Rfn=list(fn="zero"))$yp TapeR:::SK_EBLUP_LME.f(1.3/27, 30, c(1.3, 5)/27, SK.par.lme) par.lme <- SK.par.lme h <- 12 # tree height xm <- c(1.3, 3) / h # relative measuring height ym <- c(8, 7.5) # measured diameter xp <- c(0.5, 1) / h # relative prediction height TapeR:::SK_EBLUP_LME.f(xm, ym, xp, SK.par.lme)
data("SK.par.lme") TapeR:::SK_EBLUP_LME.f(1.3/27, 30, 1.3/27, SK.par.lme) ## using empirical best linear unbiased estimator: estimate != 30 TapeR:::SK_EBLUP_LME.f(1.3/27, 30, 1.3/27, SK.par.lme, Rfn=list(fn="sig2"))$yp ## interpolate / force through given diameter: estimate == 30 TapeR:::SK_EBLUP_LME.f(1.3/27, 30, 1.3/27, SK.par.lme, Rfn=list(fn="zero"))$yp TapeR:::SK_EBLUP_LME.f(1.3/27, 30, c(1.3, 5)/27, SK.par.lme) par.lme <- SK.par.lme h <- 12 # tree height xm <- c(1.3, 3) / h # relative measuring height ym <- c(8, 7.5) # measured diameter xp <- c(0.5, 1) / h # relative prediction height TapeR:::SK_EBLUP_LME.f(xm, ym, xp, SK.par.lme)
Internal function not usually called by users
SK_VOLab_EBLUP_LME.f( xm, ym, a = 0, b = 1, Ht, par.lme, Rfn = list(fn = "sig2"), IntPolOpt = TRUE, ... )
SK_VOLab_EBLUP_LME.f( xm, ym, a = 0, b = 1, Ht, par.lme, Rfn = list(fn = "sig2"), IntPolOpt = TRUE, ... )
xm |
relative heights for which measurements are available |
ym |
corresponding diameter measurements in height |
a |
relative height of lower threshold of stem section |
b |
relative height of upper threshold of stem section |
Ht |
tree height |
par.lme |
List of taper model parameters obtained by |
Rfn |
list with function name to provide estimated or assumed residual variances for the given measurements, optionally parameters for such functions |
IntPolOpt |
option for method of interpolation, if TRUE using a natural
interpolating spline ( |
... |
not currently used |
with Rfn=list(fn="zero")
one can decide whether the
measured diameters are forced to lie exactly on the taper curve; this
interferes somewhat with the IntPolOpt
, which determines the method of
taper curve point interpolation for integration. The default TRUE
(used throughout all function calls) applies natural interpolating splines,
hence this does not contradict the optional use of Rfn=list(fn="zero")
.
List with two elements, the estimated volume and its variance
Edgar Kublin
Taper model parameters for spruce in Germany based on BWI3 data obtained using
TapeR_FIT_LME.f
.
data(SK.par.lme)
data(SK.par.lme)
See Value section of TapeR_FIT_LME.f
.
Kublin, E., Breidenbach, J., Kaendler, G. (2013) A flexible stem taper and volume prediction method based on mixed-effects B-spline regression, Eur J For Res, 132:983-997.
data(SK.par.lme)
data(SK.par.lme)
Fits a taper curve model with random effects on tree-level based on B-Splines to the specified diameter-height data. Number and position of nodes and order of B-Splines can be specified.
TapeR_FIT_LME.f( Id, x, y, knt_x, ord_x, knt_z, ord_z, IdKOVb = "pdSymm", control = list(), ... )
TapeR_FIT_LME.f( Id, x, y, knt_x, ord_x, knt_z, ord_z, IdKOVb = "pdSymm", control = list(), ... )
Id |
Vector of tree identifiers of same length as diameter and height measurements. |
x |
Numeric vector of height measurements (explanatory variables) along the stem relative to the tree height. |
y |
Numeric vector of diameter measurements (response) along the stem (in centimeters). |
knt_x |
Numeric vector of relative knot positions for fixed effects. |
ord_x |
Numeric scalar. Order of fixed effects Spline (4=cubic). |
knt_z |
Numeric vector of relative knot positions for random effects. |
ord_z |
Numeric scalar. Order of random effects Spline (4=cubic). |
IdKOVb |
Character string. Type of covariance matrix used by
|
control |
a list of control values for the estimation algorithm to
replace the default values returned by the function
|
... |
not currently used |
If too few trees are given, the linear mixed model (lme) will not converge. See examples for a suggestion of node positions.
The variance parameters theta
are stored in the natural parametrization
(Pinheiro and Bates (2004), p. 93). This means log for variances and logit for
covariances. theta
is the vectorized triangle of the random effects
covariance matrix + the residual variance (lSigma). Given there are 2 inner
knots for random effects, the structure will be
c(sig^2_b1, sig_b1 sig_b2, sig_b1 sig_b3, sig_b1 sig_b4, sig^2_b2,...,sig^2_b4, lSigma)
List of model properties
fit.lmeSummary of the fitted lme model.
par.lmeList of model parameters (e.g., coefficients and
variance-covariance matrices) needed for volume estimation and other
functions in this package.
Components of the par.lme
list
knt_xRelative positions of the fixed effects Spline knots along the stem.
pad_knt_xPadded version of knt_x, as used to define B-Spline design matrix.
ord_xOrder of the spline.
knt_zRelative positions of the random effects Spline knots along the stem.
pad_knt_zPadded version of knt_z, as used to define B-Spline design matrix.
ord_zOrder of the spline.
b_fixFixed-effects spline coefficients.
KOVb_fixCovariance of fixed-effects.
sig2_epsResidual variance.
dfResResidual degrees of freedom.
KOVb_rndCovariance of random effects.
thetaVariance parameters in natural parametrization. See Details.
KOV_thetaApproximate asymptotic covariance matrix of variance parameters.
Edgar Kublin
Kublin, E., Breidenbach, J., Kaendler, G. (2013) A flexible stem taper and volume prediction method based on mixed-effects B-spline regression, Eur J For Res, 132:983-997.
E_DHx_HmDm_HT.f
, E_DHx_HmDm_HT_CIdHt.f
,
E_HDx_HmDm_HT.f
, E_VOL_AB_HmDm_HT.f
# load example data data(DxHx.df) # prepare the data (could be defined in the function directly) Id = DxHx.df[,"Id"] x = DxHx.df[,"Hx"]/DxHx.df[,"Ht"]#calculate relative heights y = DxHx.df[,"Dx"] # define the relative knot positions and order of splines knt_x = c(0.0, 0.1, 0.75, 1.0); ord_x = 4 # B-Spline knots: fix effects; order (cubic = 4) knt_z = c(0.0, 0.1 ,1.0); ord_z = 4 # B-Spline knots: rnd effects # fit the model taper.model <- TapeR_FIT_LME.f(Id, x, y, knt_x, ord_x, knt_z, ord_z, IdKOVb = "pdSymm") ## save model parameters for documentation or dissimination ## parameters can be load()-ed and used to predict the taper ## or volume using one or several measured dbh #spruce.taper.pars <- taper.model$par.lme #save(spruce.taper.pars, file="spruce.taper.pars.rdata")
# load example data data(DxHx.df) # prepare the data (could be defined in the function directly) Id = DxHx.df[,"Id"] x = DxHx.df[,"Hx"]/DxHx.df[,"Ht"]#calculate relative heights y = DxHx.df[,"Dx"] # define the relative knot positions and order of splines knt_x = c(0.0, 0.1, 0.75, 1.0); ord_x = 4 # B-Spline knots: fix effects; order (cubic = 4) knt_z = c(0.0, 0.1 ,1.0); ord_z = 4 # B-Spline knots: rnd effects # fit the model taper.model <- TapeR_FIT_LME.f(Id, x, y, knt_x, ord_x, knt_z, ord_z, IdKOVb = "pdSymm") ## save model parameters for documentation or dissimination ## parameters can be load()-ed and used to predict the taper ## or volume using one or several measured dbh #spruce.taper.pars <- taper.model$par.lme #save(spruce.taper.pars, file="spruce.taper.pars.rdata")
Internal function not usually called by users
TransKnots(knots = c(seq(0, 1, 0.1)), ord = 4, ...)
TransKnots(knots = c(seq(0, 1, 0.1)), ord = 4, ...)
knots |
knot positions for spline function |
ord |
order of the spline function |
... |
not currently used |
transformed knots vector, especially with repeated first and last knot given order of spline function
Edgar Kublin
Internal function not usually called by users
updateKOV(var, kov)
updateKOV(var, kov)
var |
vector of new variances |
kov |
variance-covariance matrix to be updated |
Firstly, kov
is transformed to a correlation matrix.
Secondly, this correlation matrix is backtransformed to a variance-covariance
matrix using the new variances.
updated variance-covariance matrix
Christian Vonderach
Internal function not usually called by users
xy0_root.f(x, y0, SK, par.lme, ...)
xy0_root.f(x, y0, SK, par.lme, ...)
x |
relative height |
y0 |
diameter for which height is required |
SK |
return of |
par.lme |
List of taper model parameters obtained by |
... |
not currently used |
used in xy0_SK_EBLUP_LME.f
to find the root of taper curve
(i.e. height x
) at given diameter y0
difference between actual diameter at height x
and given
diameter y0
Edgar Kublin
Internal function not usually called by users
xy0_SK_EBLUP_LME.f(xm, ym, y0, par.lme, Rfn = list(fn = "sig2"), ...)
xy0_SK_EBLUP_LME.f(xm, ym, y0, par.lme, Rfn = list(fn = "sig2"), ...)
xm |
relative heights for which measurements are available |
ym |
corresponding diameter measurements in height |
y0 |
given diameter for which height is required |
par.lme |
Fitted model object, return of |
Rfn |
list with function name to provide estimated or assumed residual variances for the given measurements, optionally parameters for such functions |
... |
not currently used |
function used to transform given diameter in volume calculation into
height; c.f E_VOL_AB_HmDm_HT.f
; with Rfn
one can decide
whether the measured diameters are forced to lie exactly on the taper curve
Rfn$fn="sig2"
or not Rfn$fn="zero"
. Other options are possible,
see also SK_EBLUP_LME.f
and resVar
.
relative height of given diameter y0
Edgar Kublin
Internal function not usually called by users
XZ_BSPLINE.f(x, knt, ord, ...)
XZ_BSPLINE.f(x, knt, ord, ...)
x |
relative height measurements |
knt |
knot positions for B-Splines, usually taken from a model fit by
|
ord |
order of B-Splines, usually taken from a model fit by
|
... |
not currently used |
List with height measurements (x
), the fixed effects B-splines
matrix and the random effects B-splines matrix.
Edgar Kublin
Internal function not usually called by users
y2x_isp.f(x, x.grd, y.grd, ...)
y2x_isp.f(x, x.grd, y.grd, ...)
x |
relative height |
x.grd |
relative heights for interpolation |
y.grd |
diameter of taper curve at relative heights |
... |
not currently used |
squared estimated diameter based on natural interpolating spline
(splinefun
)
Edgar Kublin
Internal function not usually called by users
y2x_ssp.f(x, x.grd, y.grd, ...)
y2x_ssp.f(x, x.grd, y.grd, ...)
x |
relative height |
x.grd |
relative heights for interpolation |
y.grd |
diameter of taper curve at relative heights |
... |
not currently used |
squared estimated diameter based on smoothing splines
(smooth.spline
)
Edgar Kublin
Internal function not usually called by users
yx_isp.f(x, x.grd, y.grd, ...)
yx_isp.f(x, x.grd, y.grd, ...)
x |
relative height |
x.grd |
relative heights for interpolation |
y.grd |
diameter of taper curve at relative heights |
... |
not currently used |
estimated diameter based on natural interpolating spline
(splinefun
)
Edgar Kublin
Internal function not usually called by users
yx_ssp.f(x, x.grd, y.grd, ...)
yx_ssp.f(x, x.grd, y.grd, ...)
x |
relative height |
x.grd |
relative heights for interpolation |
y.grd |
diameter of taper curve at relative heights |
... |
not currently used |
estimated diameter based on smoothing splines
(smooth.spline
)
Edgar Kublin