Package 'lmfor'

Title:	Functions for Forest Biometrics
Description:	Functions for different purposes related to forest biometrics, including illustrative graphics, numerical computation, modeling height-diameter relationships, prediction of tree volumes, modelling of diameter distributions and estimation off stand density using ITD. Several empirical datasets are also included.
Authors:	Lauri Mehtatalo and Kasper Kansanen
Maintainer:	Lauri Mehtatalo <[email protected]>
License:	GPL-2
Version:	1.6
Built:	2025-02-22 03:54:37 UTC
Source:	https://github.com/cran/lmfor

Help Index

Functions of Lauri Mehtatalo
Increment core data of Scots pine trees
Individual tree characteristics and ALS data
Breaking resistance (=bending strength) of birch wood samples
Plot circles of a specified radius
Evaluate the fit of a tree diameter distribution
Fit a Height-Diameter model to forest tree data using functions of package nlme.
CO2 exchange of transplanted Sphagnum fuscum moss in a chronosequence of mires.
Available 2- and 3- parameter H-D model functions to be used by function fithd.
Estimate stand density using a Horvitz–Thompson-like estimator
Estimate forest characteristics of interest in circular plot sampling using a Horvitz–Thompson-like estimator
Impute missing tree heights into a forest data using a nonlinear (mixed-effects) model.
Wood-decaying fungi carried by bark beetle individuals and their mites.
A spaghetti plot of grouped data
The Four-parameter Logit-logistic Distribution
A whiskers type residual plot
Solve a Nonlinear Equation Using Newton-Raphson algorithm.
Solve a Systems of Nonlinear Equations Using the Newton's Method
Increment core data of Scots pine trees
The Percentile-based Distribution
Sapling counts from sample plots of sapling stands in Finland.
Sapling counts from sample plots of sapling stands in Finland.
Diagnostic plot a Height-Diameter model residuals
Individual tree volume functions for Finland
Normal QQ-plot of a fitted H-D model
Properties of sample quantiles from a tree population described by the percentile-based diameter distribution.
Recovery of Weibull parameters of tree diameter distribution using measured stand characteristics
The Weibull scale parameter for the given mean/median diameter and shape parameter.
Raw sample plot data of Scots pine in Ilomantsi, Finland.
Heights and diameters of Scots pine trees in Ilomantsi, Finland.
Productivity of stump lifting machines.
Effect of thinning on individual tree growth for 62 trees.
Effect of thinning on individual tree growth
Individual tree volume modeling data for 8508 pine, spruce and birch trees in Finland
Solve a simple equation using a step halving algorithm.
The Variable Form-Factor Volume Model

Functions for different purposes related to Forest biometrics, including illustrative graphics, numerical computation, modeling height-diameter relationships, prediction of tree volumes, modeling of diameter distributions and estimation off stand density using ITD. Several empirical datasets are also included; those data sets are used in the examples of Mehtatalo and Lappi (2020a, 2020b).

Details

Package:	lmfor
Type:	Package
Version:	1.5
Date:	2020-06-16
License:	GPL -2
LazyLoad:	yes

Author(s)

Lauri Mehtatalo <[email protected]> and Kasper Kansanen <[email protected]>

References

Mehtatalo, Lauri and Lappi, Juha 2020a. Biometry for Forestry and Environmental Data: with examples in R. New York: Chapman and Hall/CRC. 426 p. doi:10.1201/9780429173462.

Mehtatalo, Lauri and Lappi, Juha 2020b. Biometry for Forestry and Environmental Data: with examples in R. Full Versions of The Web Examples. Available at http://www.biombook.org.

Increment core data of Scots pine trees

Description

Post-thinning growth ring measurements of 88 trees of a long-term thinning experiment on a naturally regenerated Scots pine stand in Eastern Finland.

Usage

data(afterthin)data(afterthin)

Format

A data frame with 1319 observations on the following 7 variables.

Plot: Sample plot id, a factor with 10 levels.
Tree: Tree id, a factor with 55 levels (same tree id may occur on different plots!).
Year: Calendar year of the ring.
SDAfterThin: Stand density (trees per ha) of the sample plot.
SDClass: Thinning treatment, factor with 4 levels (1=Control, 2=Light, 3=Moderate, 4=Heavy).
CA: Current tree age in years.
RBA: Ring Basal area, $mm^2$

Details

Long-term thinning experiment on a naturally regenerated Scots pine stand in Eastern Finland. The experiment consists of 10 sample plots, in four different classes according to the post-thinning stand density. The plots were thinned in winter 1986-1987. In winter 2006 -2007, 10 trees were felled from each plot. A radial 5mm by 5mm segment from pith to bark was cut from each tree at height 1.3 meter height. Ring widths from pith to bark were analyzed for each sample, using an ITRAX X-ray microdensitometer an post-processed to create ring widths from pith to bark were determined for each disc. The ring widths were further transformed to ring basal areas by assuming circular, growth rings. For 12 trees, ring widths could not be extracted. The data includes ring widths for a total of 88 trees between years 1991-2005. The original data is available in data set patti.

References

Mehtatalo, L., Peltola, H., Kilpelainen, A. and Ikonen, V.-P. 2014. The response of basal area growth of Scots pine to thinning: A longitudinal analysis of tree-specific series using a nonlinear mixed-effects model. Forest Science 60 (4): pp. 636-644. DOI: doi:10.5849/forsci.13-059.

Mehtatalo, Lauri and Lappi, Juha 2020. Biometry for Forestry and Environmental Data: with examples in R. New York: Chapman and Hall/CRC. 426 p. doi:10.1201/9780429173462

Examples

data(afterthin)
par(mfcol=c(2,1),cex=0.7,mai=c(0.8,0.8,0.5,0.1))
linesplot(afterthin$CA,
		  afterthin$RBA,
		  group=afterthin$Plot:afterthin$Tree,
		  col.lin=as.numeric(afterthin$SDClass),cex=0,
		  xlab="Tree age",
		  ylab=expression("Ring basal area, "*mm^2))
		  
linesplot(afterthin$Year,
		  afterthin$RBA,
		  group=afterthin$Plot:afterthin$Tree,
		  col.lin=as.numeric(afterthin$SDClass),cex=0,
		  xlab="Year",
		  ylab=expression("Ring basal area, "*mm^2))
data(afterthin)
par(mfcol=c(2,1),cex=0.7,mai=c(0.8,0.8,0.5,0.1))
linesplot(afterthin$CA,
		  afterthin$RBA,
		  group=afterthin$Plot:afterthin$Tree,
		  col.lin=as.numeric(afterthin$SDClass),cex=0,
		  xlab="Tree age",
		  ylab=expression("Ring basal area, "*mm^2))
		  
linesplot(afterthin$Year,
		  afterthin$RBA,
		  group=afterthin$Plot:afterthin$Tree,
		  col.lin=as.numeric(afterthin$SDClass),cex=0,
		  xlab="Year",
		  ylab=expression("Ring basal area, "*mm^2))

Individual tree characteristics and ALS data

Description

Field-measured and remotely sensed characteristics of 1510 individual Scots Pine trees from 56 sample plots in Kiihtelysvaara, Eastern Finland.

Usage

data(alsTree)data(alsTree)

Format

A data frame with 1510 observations on the following 15 variables.

plot: Sample plot id, integer
tree: Tree id, integer
DBH: Tree diameter at breast height, cm
H: Tree height, m
V: Tree volume, m^2)
HDB: Height of lowest dead branch, m
HCB: Crown base height, m
hmax: Maximum return height, m
h20: 20th percentile of return heights within the tree crown, m
h30: 30th percentile of return heights within the tree crown, m
h70: 70th percentile of return heights within the tree crown, m
h80: 80th percentile of return heights within the tree crown, m
a_hmean: Mean height of returns in the 250m^2 neighbourhood of the tree, m
a_veg: Proportion of returns from vegetation in the 250m^2 neighbourhood of the tree, m
a_h30: 30th percentile of returns in the 250m^2 neighbourhood around the tree, m
a_h70: 70th percentile of returns in the 250m^2 neighbourhood around the tree, m

Details

Field measurements of tree diameter and height, height of dead branch and crown base height and tree location were taken from the trees in a field campaign. Volume was estimated based on diameter, height and upper stem diameter. In addition, the area was remotely sensed using airborne laser scanning. Detectable individual trees were delineated from the ALS point cloud and associated with the field measurements. From a large set of tree-specific ALS characteristics, the data includes those that were used in the final models of stand characteristics in Maltamo et al (2012).

References

Maltamo, M., Mehtatalo, L., Vauhkonen, J. and Packalen, P. 2012. Predicting and calibrating tree attributes by means of airborne laser scanning and field measurements. Canadian Journal of Forest Research 42: 1896-1907. doi:10.1139/x2012-134

Mehtatalo, Lauri and Lappi, Juha 2020. Biometry for Forestry and Environmental Data: with examples in R. New York: Chapman and Hall/CRC. 426 p. doi:10.1201/9780429173462.

Examples

data(alsTree)
data(alsTree)

Breaking resistance (=bending strength) of birch wood samples

Description

Measurements of the breaking resistance of wood samples from a total of 118 downy birch (Betula pubescens) trees.

Usage

data(BrkRes)data(BrkRes)

Format

A data frame with 274 observations on the following 7 variables.

Tree: Tree id (numeric)
Resistance: Breaking resistance, MPa
Density: Wood density (as air-dry in 12-15% moisture) , g/cm^3
FibreLength: Fibre length, mm
RingClass: Categorical with three levels indicating the position within the stem. 1=near the pitch (inside), 2=middle, 3=near the bark (outer).
SeedOrigin: Binary variable about the origin of the tree: 1=from seed, 0=sprouted.
Site: Categorical site class with four levels.

Details

A total of 1 - 4 wood samples per tree were collected. The samples were classified to three ring classes according to the distance of the sample from the pith. Each sample was measured destructively for breaking resistance in the laboratory. The measured variables also included wood density. The data has been collected by Hanna Joronen, Katri Luostarinen and Veikko Mottonen.

References

Mehtatalo, Lauri and Lappi, Juha 2020. Biometry for Forestry and Environmental Data: with examples in R. New York: Chapman and Hall/CRC. 426 p. doi:10.1201/9780429173462

Joronen, H. 2020. Taivutusmurtolujuuteen ja kimmokertoimeen vaikuttavat tekijat hieskoivulla nuorpuussa ja aikuispuussa seka niiden mallintaminen lineaarisella sekamallilla. Master's thesis, University of Eastern Finland.

Examples

data(BrkRes)

brmod1 <- lme(Resistance~RingClass+Density, 
          random=~RingClass-1|Tree, 
          data=BrkRes)
data(BrkRes)

brmod1 <- lme(Resistance~RingClass+Density, 
          random=~RingClass-1|Tree, 
          data=BrkRes)

Plot circles of a specified radius

Description

Adds circles of radii r at coordinates specified by x and y onto an existing plot.

Usage

circle(x,y,r,border="black",lty="solid",lwd=1,fill=NULL)
circle(x,y,r,border="black",lty="solid",lwd=1,fill=NULL)

Arguments

`x`, `y`, `r`	Vectors of the x- and y- coordinates of the midpoints and the associated radii. Vectors x, y and r should be of the same length
`border`, `lty`, `lwd`	the draving color, line type and line width of the perimeter line. Use border=NA to omit the perimeter.
`fill`	The color used to fill the circles. fill=NULL does not fill at all.

Value

This function is used for its side effects on the graphical display.

Author(s)

Lauri Mehtatalo <[email protected]>

Examples

plot(0,type="n",xlim=c(-2,12),ylim=c(-2,12))
#Plot on average 7 tree crowns of Weibull-distributed radius at random locations
n<-rpois(1,7)
circle(x=runif(n,0,10),
       y=runif(n,0,10),
       r=rweibull(n,6,2))

plot(0,type="n",xlim=c(-2,12),ylim=c(-2,12))
#Plot on average 7 tree crowns of Weibull-distributed radius at random locations
n<-rpois(1,7)
circle(x=runif(n,0,10),
       y=runif(n,0,10),
       r=rweibull(n,6,2))

Evaluate the fit of a tree diameter distribution

Description

A function to compare the fit of the observed tree diameter data (d) to a specified diameter distribution (density).

Usage

ddcomp(d,density="dweibull",power=0,limits=seq(0,100),limitsd=limits,plot=FALSE,...)
ddcomp(d,density="dweibull",power=0,limits=seq(0,100),limitsd=limits,plot=FALSE,...)

Arguments

`d`	numeric vector of observed diameters
`density`	either a valid name for a probability density function in R or a vector of diameter class densities for diameter classes whose limits are given in vector limitsd
`power`	the weight used in error index. Value 2 gives BA weight, 0 (default) the unweighted
`limits`	the diameter class limits to compute the error index
`limitsd`	see the description of argument `density`
`plot`	logical. Should a graph be produced to illustrate the ecdf of d and the cdf corresponding to density
`...`	additional arguments passed to function specified by a character-type density. e.g. Weibull shape and scale of if density="dweibull"

Details

The comparison is done for mean, variance and standard deviation and shape. The shape is compared by computing the sum of absolute differences (error index) in densities for the observed data and predicted density in diameter classes specified by "limits". The error index has therefore a value between 0 (complete match) and 2 (complete mismatch). The error index is computed for the predicted density as such (ei1) and to a rescaled and switched density, which has exactly same mean and variance as the given diameter data (e12).

The error index is calculated as the sum of variable $(f_{obs}-f_{pred})x^{power}$ over the diameter classes, where $x$ is the midpoint of the diaemeter class and $f_{obs}-f_{pred}$ is the difference in predicted and observed frequency. By default, $power=0$ .

Value

A list of components

`mudif`	The difference in means
`vardif`	The difference in variances
`sddif`	The difference in standard deviations
`ei1`	the error index for original predicted distribution (see details)
`ei1`	the error index for scaled predicted distribution (see details)

Author(s)

Lauri Mehtatalo <[email protected]>

Examples

# Example 
# Observed diameters 
d<-c(18.8,24.2,18.7,13.0,18.9,22.4,17.6,22.0,18.8,22.9,
     16.7,13.7,20.6,15.1,31.8,17.2,19.6,16.8,19.3,27.4,
	 23.7,18.2,19.7,18.9,23.0,21.4,23.8,22.1,24.2,20.9)

# Weibull(5,20) distribution in 1 cm classes (class limits from 0,...,60)
f<-pweibull(1:60,5,20)-pweibull(0:59,5,20)

# compare using the classified true distribution (approximate)
ddcomp(d,density=f,limitsd=0:60,limits=0:100,plot=TRUE)

# compare b specifying a Weibull dsitribution (accurate)
ddcomp(d,density="dweibull",shape=5,scale=20,plot=TRUE)
# Example 
# Observed diameters 
d<-c(18.8,24.2,18.7,13.0,18.9,22.4,17.6,22.0,18.8,22.9,
     16.7,13.7,20.6,15.1,31.8,17.2,19.6,16.8,19.3,27.4,
	 23.7,18.2,19.7,18.9,23.0,21.4,23.8,22.1,24.2,20.9)

# Weibull(5,20) distribution in 1 cm classes (class limits from 0,...,60)
f<-pweibull(1:60,5,20)-pweibull(0:59,5,20)

# compare using the classified true distribution (approximate)
ddcomp(d,density=f,limitsd=0:60,limits=0:100,plot=TRUE)

# compare b specifying a Weibull dsitribution (accurate)
ddcomp(d,density="dweibull",shape=5,scale=20,plot=TRUE)

Fit a Height-Diameter model to forest tree data using functions of package `nlme`.

Description

Fits either linear or nonlinear Height-Diameter (H-D) model into a dataset of tree heights and diameters. Possible hierarchy of the data can be taken into account through random effects. Several commonly used nonlinear two-parameter H-D functions are available. Linear functions can be used as well.

Usage

fithd(d, h, plot=c(), modelName="naslund", nranp=2,
      random=NA, varf=0, na.omit=TRUE, start=NA, bh=1.3,
	  control = list(), SubModels=NA, vfstart=0)
fithd(d, h, plot=c(), modelName="naslund", nranp=2,
      random=NA, varf=0, na.omit=TRUE, start=NA, bh=1.3,
	  control = list(), SubModels=NA, vfstart=0)

Arguments

`d`	A numerical vector of tree diameters, usually given in cm.
`h`	A numerical vector of tree heights, usually given in meters. Should be of the same length as `d`.
`plot`	A vctor of type `numeric` or `factor`, defining the groups of the data; usually the plot indices. Should be of the same length as `d` and `h`.
`modelName`	Either (i) a character vector specifying the name of the nonlinear function or (ii) the formula specifying a linear model. In case (i) the name should be one of the functions documented on the help page of `HDmodels`. In case (ii), it should be the linear formula in the form that is entered to the function `lme`, for example `model=h~d+I(d^2)-1`
`nranp`, `random`	Parameters nranp and random specify two alternative ways to specify the random effects of the model. An easy but restricted way is to use argument `nranp`. It is an integer between 0 and the number of fixed parameters, and has the the following meaning in the case of nonlinear model: If `nranp=0`, then a model without random parameters is fitted. Results to a fixed-effects model, and argument `plot` is not used. If `nranp=1`, then parameter a of a nonlinear function or the first coefficient of the linear formula is assumed to vary among plots. If `nranp=2`, then a and b or the first two terms of the linear formula are assumed to vary among plots or If `nranp=3`, then a b, and c of a three-parameter nonlinear model or three first coefficients of a linear model are assumed to vary among plots In the case of linear model, the constant (if exists) it always counted as the first term. As an alternative to nranp, argument `random` can be used to express the random part as a nlme formula, but without specification of the grouping structure. The provided formula is passed to the lme or nlme function. Argument random is always used when provided, so nranp has effect only if `random=NA` (the default).
`varf`	Numeric with values 0, 1 or 2. If 0 or FALSE, no variance function is used. If varf=1, 2 or TRUE, then the power- type variance function var(e)=sigma^2w^(2delta) is used. where weight w is the raw diameter (when varf=1 or TRUE), or w=max(1,dsd+3) (when varf=2), where dsd=(d-D)/SDD. Here d is tree diameter, D and SDD are the mean and standard deviation of diameters on the plot in question.
`na.omit`	Should missing heights be omitted. Defaults to `TRUE`.
`start`	A vector of the starting values of the parameters of the nlme fit. If NA, then the starting values are computed using the function computing the starting values (e.g., startHDnaslund, see `HDmodels`).
`bh`	The applied breast height. Defaults to 1.3 (meters).
`control`	Parameters to control of the model fitting algorithm, see `nlmeControl` for details.
`SubModels`	Implemented only for nonlinear models. A character vector of length 2 or 3, according to the number of parameters in the model. It allows submodels for parameters a, b (and c), where the parameter is explaiend by plot-specific mean diameter ("~dmean"), plot-specific standard deviation "~dsd", or diameter standardized at plot level ("~dstd"), when the predictor is (d-D)/SDD (see teh documentation of argument varf). Defaults to NA, which corresponds to no submodels, or `submodels=c("~1","~1","1")`
`vfstart`	Starting value of the power parameter delta of the variance function. Defaults to 0.

Details

Depending on the model (nonlinear or linear, mixed-effects model or marginal), the the model is fitted using one of the following functions functions of the nlme package: nlme, lme, gls or gnls.

See available H-D functions at HDmodels. The user can define her own new functions as specified at HDmodels.

Value

An object of class hdmod, inheriting from class nlme.

Author(s)

Lauri Mehtatalo <[email protected]>

References

Mehtatalo, L., Gregoire, T.G., and de Miguel, S. Modeling Height-diameter curves for height prediction. Canadian Journal of Forest Research, 45(7): 826-837, doi:10.1139/cjfr-2015-0054

Examples

data(spati)

fithd(spati$d,spati$h,spati$plot)
fithd(spati$d,spati$h,spati$plot,SubModels=c("dmean","log(dmean)"),varf=2)
data(spati)

fithd(spati$d,spati$h,spati$plot)
fithd(spati$d,spati$h,spati$plot,SubModels=c("dmean","log(dmean)"),varf=2)

CO2 exchange of transplanted Sphagnum fuscum moss in a chronosequence of mires.

Description

The net carbon dioxide exchange of late successional moss species (Sphagnum fuscum) samples under seven levels of photosynthetic photon flux density in cronosequence of land uplift mires on the Finnish side of Bothnia Bay in Siikajoki, Finland. Moss samples were transplanted from the late succession site (Site 6) to all sites and photosynthetic activity was measured one year later for those samples which had survived.

Usage

data(foto)data(foto)

Format

A data frame with 455 observations on the following 8 variables.

Site: a factor with levels 1,...,6 from the earliest successional stage to the latest
Treatment: a factor with value 1 for samples with competitor removal treatment and 0 for untreated control.
sample: a factor with unique value for each of the 72 survived samples
moisture: moisture of the sample
PARtop: Photosyntetically active radiation (photon flux density $\mu mol /m^2 /s^2$ )
WT: water table, cm
A: Net CO2 exchange, $\mu {mol} / g / h$ ,
subplot: a factor with unique value for each replicate

Details

The number of transplanted replicates per site was 12, with two samples per replicate. One of the samples was treated with competing vegetation removal before transplanting whereas the other was left untreated. The 12 replicates per site were planted in locations with 2 to 3 different ground water table levels. A year after the transplanting, the photosynthetic activity (A) of the survived transplanted samples was recorded using seven artifically created light conditions ranging from complete darkness (PPFD=0) to extreme light conditions (PPFD=2000) using an open, fully controlled flow- through gas exchange fluorescence measurement system (GFS-3000; Walz, Effeltrich, Germany).

References

Laine, A.M., Ehonen, S., Juurola, E., Mehtatalo, L., and Tuittila, E-S. 2015. Performance of late succession species along a chronosequence: Environment does not exclude Sphagnum fuscum from the early stages of mire development. Journal of Vegetation Science 26(2): 291-301. doi:10.1111/jvs.12231

Mehtatalo, Lauri and Lappi, Juha 2020. Biometry for Forestry and Environmental Data: with examples in R. New York: Chapman and Hall/CRC. 426 p. doi:10.1201/9780429173462

Examples

data(foto)

LightResp<-function(PPFD,alpha=0.1,Pmax=10,A0=0) {
	       A0+Pmax*PPFD/(alpha+PPFD)
           }
               
library(nlme)

model5<-nlme(A~LightResp(PARtop,alpha,Pmax,A0),
		fixed=list(alpha~Site+Treatment+moisture,Pmax~Site+Treatment,A0~Site),
		random=list(sample=Pmax+alpha~1),
		data=foto,
		start=c(c(80,0,0,0,0,0,0,0),c(100,0,0,0,0,0,0),c(-20,0,0,0,0,0)),
		verbose=TRUE)
data(foto)

LightResp<-function(PPFD,alpha=0.1,Pmax=10,A0=0) {
	       A0+Pmax*PPFD/(alpha+PPFD)
           }
               
library(nlme)

model5<-nlme(A~LightResp(PARtop,alpha,Pmax,A0),
		fixed=list(alpha~Site+Treatment+moisture,Pmax~Site+Treatment,A0~Site),
		random=list(sample=Pmax+alpha~1),
		data=foto,
		start=c(c(80,0,0,0,0,0,0,0),c(100,0,0,0,0,0,0),c(-20,0,0,0,0,0)),
		verbose=TRUE)

Available 2- and 3- parameter H-D model functions to be used by function `fithd`.

Description

Nonlinear functions for modeling tree height on diameter. Usually called using fithd.

Usage

HDnaslund(d, a, b, bh=1.3) 
HDcurtis(d, a, b, bh=1.3) 
HDmichailoff(d, a, b, bh=1.3) 
HDmeyer(d, a, b, bh=1.3) 
HDpower(d, a, b, bh=1.3)
HDnaslund2(d, a, b, bh=1.3)
HDnaslund3(d, a, b, bh=1.3)
HDnaslund4(d, a, b, bh=1.3)
HDmicment(d, a, b, bh=1.3) 
HDmicment2(d, a, b, bh=1.3) 
HDwykoff(d, a, b, bh=1.3) 


HDprodan(d, a, b, c, bh=1.3) 
HDlogistic(d, a, b, c, bh=1.3) 
HDrichards(d, a, b, c, bh=1.3) 
HDweibull(d, a, b, c, bh=1.3) 
HDgomperz(d, a, b, c, bh=1.3) 
HDsibbesen(d, a, b, c, bh=1.3) 
HDkorf(d, a, b, c, bh=1.3) 
HDratkowsky(d, a, b, c, bh=1.3) 
HDhossfeldIV(d, a, b, c, bh=1.3)

startHDnaslund(d, h, bh=1.3) 
startHDcurtis(d, h, bh=1.3) 
startHDmichailoff(d, h, bh=1.3) 
startHDmeyer(d, h, bh=1.3) 
startHDpower(d, h, bh=1.3)
startHDnaslund2(d, h, bh=1.3) 
startHDnaslund3(d, h, bh=1.3) 
startHDnaslund4(d, h, bh=1.3) 
startHDmicment(d, h, bh=1.3) 
startHDmicment2(d, h, bh=1.3) 
startHDwykoff(d, h, bh=1.3) 

startHDprodan(d, h, bh=1.3) 
startHDlogistic(d, h, bh=1.3) 
startHDrichards(d, h, bh=1.3, b=0.04) 
startHDweibull(d, h, bh=1.3) 
startHDgomperz(d, h, bh=1.3) 
startHDsibbesen(d, h, bh=1.3, a=0.5) 
startHDkorf(d, h, bh=1.3) 
startHDratkowsky(d, h, bh=1.3, c=5) 
startHDhossfeldIV(d, h, bh=1.3, c=5)
HDnaslund(d, a, b, bh=1.3) 
HDcurtis(d, a, b, bh=1.3) 
HDmichailoff(d, a, b, bh=1.3) 
HDmeyer(d, a, b, bh=1.3) 
HDpower(d, a, b, bh=1.3)
HDnaslund2(d, a, b, bh=1.3)
HDnaslund3(d, a, b, bh=1.3)
HDnaslund4(d, a, b, bh=1.3)
HDmicment(d, a, b, bh=1.3) 
HDmicment2(d, a, b, bh=1.3) 
HDwykoff(d, a, b, bh=1.3) 


HDprodan(d, a, b, c, bh=1.3) 
HDlogistic(d, a, b, c, bh=1.3) 
HDrichards(d, a, b, c, bh=1.3) 
HDweibull(d, a, b, c, bh=1.3) 
HDgomperz(d, a, b, c, bh=1.3) 
HDsibbesen(d, a, b, c, bh=1.3) 
HDkorf(d, a, b, c, bh=1.3) 
HDratkowsky(d, a, b, c, bh=1.3) 
HDhossfeldIV(d, a, b, c, bh=1.3)

startHDnaslund(d, h, bh=1.3) 
startHDcurtis(d, h, bh=1.3) 
startHDmichailoff(d, h, bh=1.3) 
startHDmeyer(d, h, bh=1.3) 
startHDpower(d, h, bh=1.3)
startHDnaslund2(d, h, bh=1.3) 
startHDnaslund3(d, h, bh=1.3) 
startHDnaslund4(d, h, bh=1.3) 
startHDmicment(d, h, bh=1.3) 
startHDmicment2(d, h, bh=1.3) 
startHDwykoff(d, h, bh=1.3) 

startHDprodan(d, h, bh=1.3) 
startHDlogistic(d, h, bh=1.3) 
startHDrichards(d, h, bh=1.3, b=0.04) 
startHDweibull(d, h, bh=1.3) 
startHDgomperz(d, h, bh=1.3) 
startHDsibbesen(d, h, bh=1.3, a=0.5) 
startHDkorf(d, h, bh=1.3) 
startHDratkowsky(d, h, bh=1.3, c=5) 
startHDhossfeldIV(d, h, bh=1.3, c=5)

Arguments

`d`	A vector of tree diameters, usually in cm
`h`	A vector of tree heights, usually in m. The observed heights should be always above or equal to `bh`.
`a`, `b`, `c`	Parameters a, b (and c for 3- parameter functions) of the applied function. See details for expressions of different functions.
`bh`	The applied height for the measurement of tree diameter (so called breast height). Of the same unit as `h`.

Details

The available 2- parameter functions are

Naslund: $h(d) = bh + \frac{d^2}{(a + bd)^2}$
Curtis: $h(d) = bh + a \left(\frac{d}{1 + d}\right)^b$
Michailoff: $h(d) = bh + a e^{-b d^{-1}}$
Meyer: $h(d) = bh + a (1-e^{-b d})$
Power: $h(d) = bh + a d^b$
Naslund2: $h(d) = bh + \frac{d^2}{\left(a + e^b d\right)^2}$
Naslund3: $h(d) = bh + \frac{d^2}{(e^a + b d)^2}$
Naslund4: $h(d) = bh + \frac{d^2}{(e^a + e^b d)^2}$
Michaelis-Menten: $h(d) = bh + \frac{a d}{b + d}$
Michaelis-Menten2: $h(d) = bh + \frac{d}{a + b * d}$
Wykoff: $h(d) = bh + \exp\left(a + \frac{b}{d + 1}\right)$

The available 3- parameter functions are

Prodan: $h(d) = bh + \frac{d^2}{a + bd + c d^2}$
Logistic: $h(d) = bh + \frac{a}{1 + b e^{-c d}}$
Chapman-Richards: $h(d) = bh + a (1 - e^{-bd})^c$
Weibull: $h(d) = bh + a (1 - e^{-b d^c})$
Gomperz: $h(d) = bh + a \exp(-b \exp(-c d))$
Sibbesen: $h(d) = bh + a d^{b d^{-c}}$
Korf: $h(d) = bh + a \exp(-b d^{-c})$
Ratkowsky: $h(d) = bh + a \exp\left(\frac{-b}{d + c}\right)$
Hossfeld IV: $h(d) = bh + \frac{a}{1 + \frac{1}{bd^c}}$

For each model, two functions are provided: one computing the value of the H-D model for given diameters using given values of parameters a, b (and c), and another returning the initial guesses of a, b (and c) for given h-d data.

The initial guesses are in most cases computed by fitting a linearized version of the model into the provided h-d data using lm. For some 3- parameter versions, no straightforward linearization is possible and one of the parameters is set to a fixed sensible constant. Those values can be seen as additional arguments in the corresponding startHD - functions. Details can be seen directly from the function definitions.

The user can define her own functions to be used with fithd. The case-sensitive naming of the functions should follow exactly the naming convention shown above. In addition, the names of the of arguments, as well as their order, should be the same as in the functions above.

The models are named according to references in

Zeide, B. 1993. Analysis of growth equations. Forest Science 39(3):594-616. doi:10.1093/forestscience/39.3.594
Huang, S., Titus, S.J., and Wiens, D.P. 1992. Comparison of nonlinear height-diameter functiond for major Alberta tree species. Can J. For. Res. 22: 1297-1304. doi:10.1139/x92-172

Suggestions on naming and references on the functions are welcome.

Value

For functions HDxxx, a vector of tree heights corresponding diameters d is returned. For functions startHDxxx, a named vector of initial estimates of a, b and (c).

Author(s)

Lauri Mehtatalo <[email protected]>

References

Mehtatalo, L., Gregoire, T.G., and de Miguel, S. Modeling Height-diameter curves for height prediction. Canadian Journal of Forest Research, 45(7): 826-837, doi:10.1139/cjfr-2015-0054

Examples

data(spati)
theta<-startHDnaslund(spati$d,spati$h)
plot(spati$d,spati$h)
d<-seq(0,50)
lines(d,HDnaslund(d,a=theta[1],b=theta[2]),col="red",lwd=5)
data(spati)
theta<-startHDnaslund(spati$d,spati$h)
plot(spati$d,spati$h)
d<-seq(0,50)
lines(d,HDnaslund(d,a=theta[1],b=theta[2]),col="red",lwd=5)

Estimate stand density using a Horvitz–Thompson-like estimator

Description

HTest calculates the Horvitz–Thompson-like stand density estimate (number of trees) in a specified area based on a collection of detected trees.

area_esh is an internal function for surface area calculations that can handle empty sets.

gg_wind is an internal function that forms a union of discs based on their center points and radii.

Usage

HTest(treelist, plotwindow, alpha)

area_esh(W)

gg_wind(treelist)
HTest(treelist, plotwindow, alpha)

area_esh(W)

gg_wind(treelist)

Arguments

`treelist`	A 3-column matrix containing the x and y coordinates of detected trees and their crown radii.
`plotwindow`	A `spatstat` object of class ”`owin`”, representing the area where stand density estimation is done.
`alpha`	A tuning parameter that controls the calculation of detection probabilities, or detectabilities. Must have a value from -1 to 1.
`W`	A `spatstat` object of class ”`owin`” or NULL.

Details

HTest is the Horvitz–Thompson-like stand density estimator presented by Kansanen et al. (2016) to adjust individually detected trees for non-detection. It uses individual tree detection data, namely the locations and crown radii of detected trees, to calculate detection probabilities, or detectabilities, for every detected tree, and produces an estimate based on the detectabilities. The detectability for a certain tree is based on the planar set formed by the larger trees. The parameter alpha controls how easy it is to detect a tree of certain size from under the larger trees. If alpha=1, then the tree will be detected if it is not fully covered by the larger crowns. If alpha=0, the tree will be detected if its center point is not covered. If alpha=-1, the tree will be detected if it is fully outside the larger tree crowns.

The object treelist can include trees that are not in the estimation area specified by plotwindow. This can be useful to take into account possible edge effects, by including trees with center points outside plotwindow that have crown discs that intersect plotwindow. The estimate is calculated only using those trees that have crown center points in plotwindow.

area_esh and gg_wind are internal helper functions used by HTest. First one is a shell for the spatstat function area.owin that takes into account that an intersection of two sets can be empty, represented in the calculations as NULL. The function returns 0 in this case. Otherwise, it returns the surface area of the window W. The latter function forms a union of discs that is needed in the detectability calculations.

Value

HTest returns a list with two components:

`N`	The estimated number of trees in `plotwindow`
`treelist`	matrix with columns ”r” and ”detectability”, giving the tree crown radii that have been used in the estimation, as well as the detectabilities for trees with those crown radii.

area_esh returns 0, if W is NULL; otherwise, the surface area of W.

gg_wind returns a spatstat object of class ”owin” representing a set formed as a union of discs.

Note

These functions require the package spatstat (Baddeley et al. 2015) to work.

Author(s)

Kasper Kansanen <[email protected]>

References

Kansanen, K., Vauhkonen, J., Lahivaara, T., and Mehtatalo., L. (2016) Stand density estimators based on individual tree detection and stochastic geometry. Canadian Journal of Forest Research 46(11):1359–1366. doi:10.1139/cjfr-2016-0181.

Baddeley, A., Rubak, E. and Turner, R. (2015) Spatial Point Patterns: Methodology and Applications with R. Chapman and Hall/CRC Press, London. doi:10.1201/b19708

Kansanen, K., Packalen, P., Lahivaara, T., Seppanen, A., Vauhkonen, J., Maltamo, M., and Mehtatalo., L. (2019) Horvitz–Thompson-like stand density estimation and functional k-NN in individual tree detection. Submitted manuscript.

Examples

# Generate a 10x10 meter square window:
w<-square(10)

# Generate 6 detected trees, 5 located in the window:
x<-cbind(c(6.75, 8.65, 3.95, 2, 2, 11), 
         c(1.36, 3.10, 6.66, 2, 4, 11), 
         c(1.29, 2.31, 1.80, 2, 1.5, 3))

# Draw the set formed by  the detected tree crowns:
plot(w)
plot(gg_wind(x), add=TRUE)

# Calculate the results with different alpha:
HTest(x, w, 1)
HTest(x, w, 0)
HTest(x, w, -0.75)
# Generate a 10x10 meter square window:
w<-square(10)

# Generate 6 detected trees, 5 located in the window:
x<-cbind(c(6.75, 8.65, 3.95, 2, 2, 11), 
         c(1.36, 3.10, 6.66, 2, 4, 11), 
         c(1.29, 2.31, 1.80, 2, 1.5, 3))

# Draw the set formed by  the detected tree crowns:
plot(w)
plot(gg_wind(x), add=TRUE)

# Calculate the results with different alpha:
HTest(x, w, 1)
HTest(x, w, 0)
HTest(x, w, -0.75)

Estimate forest characteristics of interest in circular plot sampling using a Horvitz–Thompson-like estimator

Description

HTest_cps calculates Horvitz–Thompson-like estimates of forest characteristics of interest in a specified circular area based on a collection of detected trees and their detection probabilities, or detectabilities. Also produces estimated variances and confidence intervals.

detectability_cps calculates detectabilities of trees in a circular plot sample.

visibility_thinning_cps takes a tree list and determines if the trees can be detected when a certain visibility-based detection condition is used.

ordering_cps is a helper function for preprocessing of tree lists: it takes a tree list and orders the trees based on their distance to plot centre point.

polar_to_cart and cart_to_polar are internal functions for transforming polar coordinates to cartesian coordinates and vice versa.

triangle_coords is an internal function that, given locations and diameters of discs, returns coordinates needed to define the areas behind the discs that are non-visible from the origin.

shades is an internal function that forms polygonal approximations of the planar sets that are non-visible from the origin.

Usage

HTest_cps(data, total=TRUE, confidence.level=0.95)
detectability_cps(data, plot.radius, alpha=0, polar=TRUE, npoly=1024, delta=NULL)
visibility_thinning_cps(data, plot.radius, alpha=0, polar=TRUE, npoly=1024, delta=NULL)
ordering_cps(data, polar=TRUE)
polar_to_cart(X)
cart_to_polar(X)
triangle_coords(X, plot.radius, polar=TRUE)
shades(X, plot.radius, polar=TRUE)
HTest_cps(data, total=TRUE, confidence.level=0.95)
detectability_cps(data, plot.radius, alpha=0, polar=TRUE, npoly=1024, delta=NULL)
visibility_thinning_cps(data, plot.radius, alpha=0, polar=TRUE, npoly=1024, delta=NULL)
ordering_cps(data, polar=TRUE)
polar_to_cart(X)
cart_to_polar(X)
triangle_coords(X, plot.radius, polar=TRUE)
shades(X, plot.radius, polar=TRUE)

Arguments

`data`	For `HTest_cps` a matrix where each row corresponds to a tree, the first column contains the detectabilities of the trees calculated with `detectability_cps`, and other columns correspond to measurements over which estimation is wanted. For `detectability_cps` a four-column matrix, each row containing the coordinates of a tree, the diameter of the tree, and an indicator if the tree has been detected (1) or not (0). In other words, each row is a vector of the form (r, phi, diameter, detected) if the locations are given in polar coordinates, or (x, y, diameter, detected) if the locations are given in cartesian coordinates. It is assumed that the trees are a sample from circular plot sampling, and that the centre point of the plot is the origin (0,0). For `visibility_thinning_cps` a three-column matrix, otherwise similar to the matrix needed by `detectability_cps` but without the ”detected” column. For `ordering_cps` a matrix, each row corresponding to a tree, the first three columns being either (x, y, diameter) or (r, phi, diameter).
`total`	Do you want to estimate population totals (TRUE) or population means (FALSE)?
`confidence.level`	The level of the approximate confidence interval, a value between 0 and 1. For example, `confidence.level`=0.95 indicates a 95 per cent confidence interval.
`plot.radius`	Radius of the plot in which circular plot sampling has been performed.
`alpha`	A tuning parameter that controls the calculation of detection probabilities, or detectabilities.
`polar`	Are the locations of trees given in polar coordinates (TRUE) or cartesian coordinates (FALSE)?
`npoly`	Number of edges for the polygonal approximation of the plot boundary and the circles used to calculate the detection probabilities. Used if `delta`is NULL.
`delta`	The tolerance of the polygonal approximation of of the plot boundary and the circles used to calculate the detection probabilities: the length of the arc that will be replaced by one edge of the polygon. If given value that is different from NULL tihs will override `npoly`
`X`	For `polar_to_cart` a two-column matrix, each row containing the polar coordinates (r, phi) of a point. For `cart_to_polar` a two-column matrix, each row containing the cartesian coordinates (x, y) of a point. For `triangle_coords` and `shades` a three-column matrix, each row containing either the polar coordinates and a diameter of a tree (r, phi, diameter) or the cartesian coordinates and a diameter of a tree (x, y, diameter). The coordinates define the centre point of the tree.

Details

The function HTest_cps produces estimates of forest characteristics of interest in a circular plot sampling situation. More specifically, it is assumed that an observer stands in a point and observes such trees that are within the fixed-area plot and are not hiding behind other trees. It is assumed that observer can record the locations and diameters of the trees that they observe. The observer can be a person or a piece of equipment, such as terrestrial laser scanner or camera.

The estimation is based on a Horvitz–Thompson-like estimator presented by Kansanen et al. (2019). This construction uses approximated detection probabilities, or detectabilities, that depend on the size and distance from the plot centre of the tree for which the probability is calculated, the nonvisible area produced by trees that are closer to the centre point, and a visibility based detection condition. It is assumed that the centre point of the sampling plot is the origin of the plane (0,0). The function detectability_cps is used to calculate the detectabilities.

The confidence intervals that HTest_cps produces are based on the t-distribution if less than 50 trees have been observed, and the standard normal distribution otherwise.

The parameter alpha is a value between -1 and 1 and it controls the detection condition. alpha=1 means that trees are detected if the stems are fully visible to the observer, alpha=0 means that they are detected if the center point is visible, and alpha=-1 means that a tree is detected if any part of the stem is visible.

The estimation is not possible if the data contains trees that cover the plot centre point.

All of the variables related to distance and size, meaning the cartesian coordinates, the distance coordinate, plot.radius, delta, and tree diameters, should have the same unit, e.g. they should all be in metres.

The function visibility_thinning_cps is useful for simulation-based testing of the estimator. Given a tree list, it classifies trees as either detected or not detected based on a visibility based detection condition.

The function ordering_cps is a useful preprocessing step for tree lists over which estimation is needed. It reorders the rows of the tree list, corresponding to trees, based on the closest distance from the stem disc to the plot centre point. detectability_cps and visibility_thinning_cps use this ordering for their advantage, as this is the assumed order or sequence of detection. Be warned that even if you do not order the data with ordering_cps these functions will, and will output tree lists with this ordering!

The functions polar_to_cart, cart_to_polar, triangle_coords, and shades are internal functions used by the two main functions. They can be useful for visualizing data.

Value

HTest_cps returns a four-column matrix, the columns containing an estimate, estimated variance of the estimator, and lower and upper bound of the approximate confidence interval for the estimate. Rows of the matrix correspond to the forest characteristics of interest given in columns 2:ncol(data) of the input matrix data. If the input matrix has named columns, these names are used as row names of the output matrix.

detectability_cps returns a matrix with the locations and diameters of the trees given as input, the indicators of their detection, and the estimated detection probabilities.

visibility_thinning_cps returns a four-column matrix with the locations and diameters of the trees given as input and the indicators of their detection: 1, if a tree has been detected based on the detection condition given by alpha, and 0, otherwise.

ordering_cps returns a matrix with same dimensions as the input matrix, rows being reordered in the assumed order of detection.

polar_to_cart returns a two-column matrix of cartesian coordinates (x, y).

cart_to_polar returns a two-column matrix of polar coordinates (r, phi).

triangle_coords returns an eight-column matrix, each row containing cartesian coordinates needed for forming a polygonal representation of an area behind a tree, nonvisible from the origin.

shades returns a list of owin objects, each representing an area behind a tree, nonvisible from the origin.

Note

These functions require the package spatstat (Baddeley et al. 2015) to work.

Author(s)

Kasper Kansanen <[email protected]>

References

Baddeley, A., Rubak, E. and Turner, R. (2015) Spatial Point Patterns: Methodology and Applications with R. Chapman and Hall/CRC Press, London. doi:10.1201/b19708

Kansanen, K., Packalen, P., Maltamo, M., and Mehtatalo, L. (2020+) Horvitz–Thompson-like estimation with distance-based detection probabilities for circular plot sampling of forests. Biometrics. doi:10.1111/biom.13312

Examples

## Not run: 
# Simulate a plot of radius 10 metres and stem density of 1000 trees/ha from the Poisson process:
set.seed(1)
N <- rpois(1, lambda=0.1*pi*10^2)
proc <- cbind(10*sqrt(runif(N, 0 ,1)),  runif(N, -pi, pi), rweibull(N, shape=12, scale=22)/100)

# Preprocess the data to the right order:
proc <- ordering_cps(proc)

# Thin the data:
thinned<-visibility_thinning_cps(data=proc, plot.radius=10, alpha=1)

# Calculate the detection probabilities:
detprob <- detectability_cps(data=thinned, plot.radius=10, alpha=1)

# Calculate estimates of stand density (number of trees) and basal area 
# (the sum of areas covered by tree stems at breast height):
mydata<-cbind(detprob[,5], 1, pi*detprob[,3]^2)

HTest_cps(data=mydata)

# Calculate estimate of mean DBH and a 99 per cent approximate confidence interval:
mydata<-cbind(detprob[,5], detprob[,3])

HTest_cps(data=mydata, total=FALSE, confidence.level=0.99)

## End(Not run)
## Not run: 
# Simulate a plot of radius 10 metres and stem density of 1000 trees/ha from the Poisson process:
set.seed(1)
N <- rpois(1, lambda=0.1*pi*10^2)
proc <- cbind(10*sqrt(runif(N, 0 ,1)),  runif(N, -pi, pi), rweibull(N, shape=12, scale=22)/100)

# Preprocess the data to the right order:
proc <- ordering_cps(proc)

# Thin the data:
thinned<-visibility_thinning_cps(data=proc, plot.radius=10, alpha=1)

# Calculate the detection probabilities:
detprob <- detectability_cps(data=thinned, plot.radius=10, alpha=1)

# Calculate estimates of stand density (number of trees) and basal area 
# (the sum of areas covered by tree stems at breast height):
mydata<-cbind(detprob[,5], 1, pi*detprob[,3]^2)

HTest_cps(data=mydata)

# Calculate estimate of mean DBH and a 99 per cent approximate confidence interval:
mydata<-cbind(detprob[,5], detprob[,3])

HTest_cps(data=mydata, total=FALSE, confidence.level=0.99)

## End(Not run)

Impute missing tree heights into a forest data using a nonlinear (mixed-effects) model.

Description

A function to impute tree heights in a forest inventory situation where all trees have been measured for diameter but only some trees have been measured for height.

Usage

ImputeHeights(d, h, plot, modelName = "naslund", nranp = 2, varf = TRUE, 
              addResidual = FALSE, makeplot=TRUE, level = 1,
              start=NA, bh=1.3, control=list(),random=NA)
ImputeHeights(d, h, plot, modelName = "naslund", nranp = 2, varf = TRUE, 
              addResidual = FALSE, makeplot=TRUE, level = 1,
              start=NA, bh=1.3, control=list(),random=NA)

Arguments

`d`	A numerical vector of tree diameters, usually given in cm.
`h`	A numerical vector of tree heights, usually given in meters. Should be of the same length as `d`.
`plot`	A vector of type `numeric` or `factor`, defining the groups of the data; usually the plot indices. Should be of the same length as `d` and `h`.
`modelName`	Either (i) a character vector specifying the name of the nonlinear function or (ii) the formula specifying a linear model. In case (i) the name should be one of the functions documented on the help page of `HDmodels`. In case (ii), it should be the linear formula in the form that is entered to the function `lme`, for example `model=h~d+I(d^2)-1`
`nranp`	Parameters nranp and random specify two alternative ways to specify the random effects of the model. An easy but restricted way is to use argument `nranp`. It is an integer between 0 and the number of fixed parameters, and has the the following meaning in the case of nonlinear model: If `nranp=0`, then a model without random parameters is fitted. Results to a fixed-effects model, and argument `plot` is not used. If `nranp=1`, then parameter a of a nonlinear function or the first coefficient of the linear formula is assumed to vary among plots. If `nranp=2`, then a and b or the first two terms of the linear formula are assumed to vary among plots or If `nranp=3`, then a b, and c of a three-parameter nonlinear model or three first coefficients of a linear model are assumed to vary among plots In the case of linear model, the constant (if exists) it always counted as the first term. As an alternative to nranp, argument `random` can be used to express the random part as a nlme formula, but without specification of the grouping structure. The provided formula is passed to the lme or nlme function. Argument random is always used when provided, so nranp has effect only if `random=NA` (the default).
`varf`	Numeric with values 0, 1 or 2. If 0 or FALSE, no variance function is used. If varf=1, 2 or TRUE, then the power- type variance function var(e)=sigma^2w^(2delta) is used. where weight w is the raw diameter (when varf=1 or TRUE), or w=max(1,dsd+3) (when varf=2), where dsd=(d-D)/SDD. Here d is tree diameter, D and SDD are the mean and standard deviation of diameters on the plot in question.
`addResidual`	Boolean. If `TRUE`, a random residual is added to the imputed height from a normal distribution using the estimated variance function. If also `level=0` or if the plot did not include any measured heights to predict the random effects, then also a randomly selected plot effect from among the predicted plot effects is added. The added plot effect is the same for all trees of a given plot.
`makeplot`	Should a residual plot of the fitted model be produced for evaluation of goodness of fit? The plot is produced using the default arguments of function `plot.hdmod`, and is not affected by the value of arguments `level` and `addResidual`.
`level`	The level of prediction. 0 means fixed-effect prediction and 1 means plot-level prediction using the random effects. Has no effect if `nranp=0`.
`start`, `bh`, `control`, `random`	Arguments passed to `fithd`. See documentation of `fithd`.

Details

The function predicts the missing heights using a nonlinear mixed-effects model or a nonlinear fixed-effects model. In mixed-effects model, plot-specific random effects can be used if other tree heights have been measured from the same plot. Also random, normally distributed residual can be added to the heights according to the estimated constant or heteroscedastic residual variance structure.

Value

A list of components

`h`	A vector of tree heights, including the measured heights for the trees with known height and imputed heights for the others.
`imputed`	A booelan vector of the same length as h, having value TRUE for imputed heights. Produced as `is.na(data$h)`
`model`	The fitted model that was used in imputation. Fitted using `fithd` which in turn calls `nlme` or `gnls`
`predType`	A vector of the same length as h, including information on the level of prediction. Value 0 means a measured height (no model prediction is used), value 1 means the plot-level prediction has been done using the estimated plot effects. Value 2 means that no sample trees were available and the prediction is based on fixed part only (if level=0) or on a simulated plot effect (if level=1).
`hpred`	Predicted heights for all trees. Equals to vector h for trees that had missing heights.

Note

Works only with the nonlinear functions specified in HDmodels; does not work if the modelName is specified as a linear expression.

Author(s)

Lauri Mehtatalo <[email protected]>

References

Mehtatalo, L., Gregoire, T.G., and de Miguel, S. Modeling Height-diameter curves for height prediction. Canadian Journal of Forest Research, 45(7): 826-837, doi:10.1139/cjfr-2015-0054

Examples

data(spati) 

ImpFixed<-ImputeHeights(spati$d,spati$h,spati$plot,level=0)
ImpRandom<-ImputeHeights(spati$d,spati$h,spati$plot,level=1,makeplot=FALSE)
# Try also
# ImpRanRes<-ImputeHeights(spati$d,spati$h,spati$plot,level=1,addResidual=TRUE,makeplot=FALSE)

plot(spati$d[!is.na(spati$h)],
     spati$h[!is.na(spati$h)],
     col=spati$plot[!is.na(spati$h)],
     main="Observations", xlab="d, cm", ylab="h, m",
     ylim=c(0,30))

plot(spati$d[ImpFixed$imputed],
     ImpFixed$h[ImpFixed$imputed],
     col=spati$plot[ImpFixed$imputed],
     main="Imputed, Naslund, Fixed", xlab="d, cm", ylab="h, m",
     ylim=c(0,30))
       
plot(spati$d[ImpRandom$imputed],
     ImpRandom$h[ImpRandom$imputed],
     col=spati$plot[ImpRandom$imputed],
     main="Imputed, Naslund, Fixed + Plot", xlab="d, cm", ylab="h, m",
     ylim=c(0,30))

# Try also
# plot(spati$d[ImpRanRes$imputed],
#     ImpRanRes$h[ImpRanRes$imputed],
#     col=spati$plot[ImpRanRes$imputed],
#     main="Imputed, Naslund, Fixed + Plot + Tree", xlab="d, cm", ylab="h, m",
#     ylim=c(0,30))
data(spati) 

ImpFixed<-ImputeHeights(spati$d,spati$h,spati$plot,level=0)
ImpRandom<-ImputeHeights(spati$d,spati$h,spati$plot,level=1,makeplot=FALSE)
# Try also
# ImpRanRes<-ImputeHeights(spati$d,spati$h,spati$plot,level=1,addResidual=TRUE,makeplot=FALSE)

plot(spati$d[!is.na(spati$h)],
     spati$h[!is.na(spati$h)],
     col=spati$plot[!is.na(spati$h)],
     main="Observations", xlab="d, cm", ylab="h, m",
     ylim=c(0,30))

plot(spati$d[ImpFixed$imputed],
     ImpFixed$h[ImpFixed$imputed],
     col=spati$plot[ImpFixed$imputed],
     main="Imputed, Naslund, Fixed", xlab="d, cm", ylab="h, m",
     ylim=c(0,30))
       
plot(spati$d[ImpRandom$imputed],
     ImpRandom$h[ImpRandom$imputed],
     col=spati$plot[ImpRandom$imputed],
     main="Imputed, Naslund, Fixed + Plot", xlab="d, cm", ylab="h, m",
     ylim=c(0,30))

# Try also
# plot(spati$d[ImpRanRes$imputed],
#     ImpRanRes$h[ImpRanRes$imputed],
#     col=spati$plot[ImpRanRes$imputed],
#     main="Imputed, Naslund, Fixed + Plot + Tree", xlab="d, cm", ylab="h, m",
#     ylim=c(0,30))

Wood-decaying fungi carried by bark beetle individuals and their mites.

Description

The total number of fungal species (ophistomatoid and non-ophistomatoid fungi are coded separately) associated with Ips typographus bark beetle individuals and their mites.

Usage

data(ips)data(ips)

Format

A data frame with 298 observations (bark betle individuals) on the following 5 variables.

Fungi: The total number of fungal species associated with the individual bark beetle.
Ophi: The number of ophistomatoid fungal species.
Other: The number of non-ophistomatoid fungal species. The three first variables are related through $Other+Ophi=Fungi$ .
Season: Categorical time of data collection with three levels: spring, summer or fall. The default is spring.
Mites: The number of mites found in the bark beetle.

Details

The ophiostomatoid fungal families Microascales and Ophiostomatales are common associates of bark beetle Ips typographus, which they use to spread within the wood material. The number of fungal species in these families is high, and a certain beetle individual can carry several fungal species with it. The bark beetles may have mites attached to them, and it may be possible that some fungal species are associated to the beetles only through the mites.

The dataset includes measurements of 289 bark beetle individuals from a storm-felled Norway spruce forest in eastern Finland. For each individual, the number of attached mites was determined using a microscope. In addition the number of fungal species per bark beetle was determined genetically. However, it was not possible to determine whether the fungi were associated with the mites or the bark beetle itself. The observations were collected at three different seasons: spring, summer and fall of the same year, approximately 100 individuals in each season. The data are used to analyze the effects of season and number of mites on the number of fungal species per bark beetle.

References

Linnakoski, R., Mahilainen, S., Harrington, A., Vanhanen, H., Eriksson, M., Mehtatalo, L., Pappinen, A., Wingfield, M.J. 2016. The seasonal succession of fungi associated with Ips typographus beetles and their phoretic mites in an outbreak region of Finland. PLOS ONE. doi:10.1371/journal.pone.0155622.

Mehtatalo, Lauri and Lappi, Juha 2020. Biometry for Forestry and Environmental Data: with examples in R. New York: Chapman and Hall/CRC. 426 p. doi:10.1201/9780429173462

Examples

data(ips)

ips$Mites2<-ips$Mites-mean(ips$Mites)

mod1<-glm(Fungi~Season+Mites,family=poisson,data=ips)
data(ips)

ips$Mites2<-ips$Mites-mean(ips$Mites)

mod1<-glm(Fungi~Season+Mites,family=poisson,data=ips)

A spaghetti plot of grouped data

Description

Orders the observations by x and thereafter plots y on x and connects observations of the same group by lines. Useful, for example, to plot a longitudinal dataset.

Usage

linesplot(x, y, group, xlab = "x", ylab = "y", 
          main = "", cex = 0.5, pch = 19, col = 1, col.lin = 1, 
          lw = FALSE, ylim = NULL, xlim = NULL, add = FALSE, lty = "solid", lwd=1)
linesplot(x, y, group, xlab = "x", ylab = "y", 
          main = "", cex = 0.5, pch = 19, col = 1, col.lin = 1, 
          lw = FALSE, ylim = NULL, xlim = NULL, add = FALSE, lty = "solid", lwd=1)

Arguments

`x`, `y`	Numerical vectors of the same length including the x and y variables.
`group`	The variable specifying the group. Should be of the same length as vectors x and y.
`xlab`, `ylab`, `main`, `cex`, `pch`, `col`, `col.lin`, `xlim`, `ylim`, `lty`, `lwd`	Graphical parameters, see `par`
`lw`	Boolean. Whether a loess smoother to be added onto the plot.
`add`	Boolean. Whether to add to an existing plot or to open a new window.

Details

The observations within the group are connected at the increasing order of x.

Value

Used for its side effects.

Author(s)

Lauri Mehtatalo <[email protected]>

References

Mehtatalo, Lauri and Lappi, Juha 2020. Biometry for Forestry and Environmental Data: with examples in R. New York: Chapman and Hall/CRC. 426 p. doi:10.1201/9780429173462

Examples

D<-rep(seq(10,30),10)
H<-(20+rep(rnorm(10,0,0.5),each=21))*exp(-1.5*D^(-1.3))
plot<-rep(1:10,each=21)
linesplot(D,H,plot)
D<-rep(seq(10,30),10)
H<-(20+rep(rnorm(10,0,0.5),each=21))*exp(-1.5*D^(-1.3))
plot<-rep(1:10,each=21)
linesplot(D,H,plot)

The Four-parameter Logit-logistic Distribution

Description

Density, distribution function, quantile function and random generation for the four-parameter logit-logistic distribution.

Usage

dll(x, mu, sigma, xi=0, lambda=1, log = FALSE)
pll(q, mu, sigma, xi=0, lambda=1, lower.tail=TRUE, log.p=FALSE)
qll(p, mu, sigma, xi=0, lambda=1, lower.tail=TRUE, log.p=FALSE)
rll(n, mu, sigma, xi=0, lambda=1)
dll(x, mu, sigma, xi=0, lambda=1, log = FALSE)
pll(q, mu, sigma, xi=0, lambda=1, lower.tail=TRUE, log.p=FALSE)
qll(p, mu, sigma, xi=0, lambda=1, lower.tail=TRUE, log.p=FALSE)
rll(n, mu, sigma, xi=0, lambda=1)

Arguments

`x`, `q`	vector of quantiles
`p`	vector of probabilitiies
`n`	number of observations. If `length(n) > 1`, the length is taken to be the number required.
`mu`, `sigma`, `xi`, `lambda`	parameters of the distribution, `xi` (minimum) defaults to 0 and `lambda` (max-min) to 1.
`log`, `log.p`	logical; if TRUE, probabilities p are given as log(p)
`lower.tail`	logical; if TRUE (default), probabilities are $P[X \le x]$ otherwise $P[X > x]$ .

Details

The logit-logistic cdf and pdf are

$F(d|\xi,\lambda,\mu,\sigma) = \frac{1}{1+e^{(\frac{\mu}{\sigma})} (\frac{d-\xi}{\xi+\lambda-d})^{-\frac{1}{\sigma}}}$

$f(d|\xi,\lambda,\mu,\sigma) = \frac{\lambda}{\sigma}\frac{1}{(d-\xi)(\xi+\lambda-d)} \frac{1}{e^{-\frac{\mu}{\sigma}}(\frac{d-\xi}{\xi+\lambda-d})^{\frac{1}{\sigma}}+e^{\frac{\mu}{\sigma}}(\frac{d-\xi}{\xi+\lambda-d})^{-\frac{1}{\sigma}}+2}$

Parameter $\xi$ is the minimum, $\lambda>0$ the width of range (max-min), $\mu$ controls the skewness and $\sigma$ the curtosis.

Value

dll gives the density, pll gives the distribution function, qll gives the quantile function, and rll generates random deviates.

Invalid arguments will result in return value NaN.

The length of the result is determined by n for rll, and is the maximum of the lengths of the numerical arguments for the other functions.

The numerical arguments other than n are recycled to the length of the result. Only the first elements of the logical arguments are used.

Author(s)

Lauri Mehtatalo <[email protected]>

References

Mingliang Wang and Keith Rennolls, 2005. Tree diameter distribution modelling: introducing the logit-logistic distribution. Canadian Journal of Forest Research, 35(6): 1305-1313, doi:10.1139/x05-057.

Mehtatalo, Lauri and Lappi, Juha 2020. Biometry for Forestry and Environmental Data: with examples in R. New York: Chapman and Hall/CRC. 426 p. doi:10.1201/9780429173462

Examples

data(spati)
d<-spati$d[spati$plot==22]
hist(d,freq=FALSE)
d0<-seq(0,60,0.1)
lines(d0,dll(d0,0.630,0.573,3.561,35.2))
data(spati)
d<-spati$d[spati$plot==22]
hist(d,freq=FALSE)
d0<-seq(0,60,0.1)
lines(d0,dll(d0,0.630,0.573,3.561,35.2))

A whiskers type residual plot

Description

A function for adding vertical lines onto residual plots to show

95% confidence intervals of means or
95% confidence intervals for individual observations

in the classes of the variable on the x-axis. Plot of the first type is useful for analyzing the fit of the assumed fixed part and plots of type b can be used to analyze the homogeneity of residuals.

Usage

mywhiskers(x, y, 
           nclass = 10, 
           limits = NA, 
           add = FALSE, 
           se = TRUE, 
           main = "", 
           xlab = "x", 
           ylab = "y", 
           ylim = NA, 
           lwd = 1, 
           highlight = "red")
mywhiskers(x, y, 
           nclass = 10, 
           limits = NA, 
           add = FALSE, 
           se = TRUE, 
           main = "", 
           xlab = "x", 
           ylab = "y", 
           ylim = NA, 
           lwd = 1, 
           highlight = "red")

Arguments

`x`	The variable on the x-axis. Usually one of the predictors or the predicted value.
`y`	The variable on the y-axis. Usually model residual.
`nclass`	The maximum number of classes to be used.
`limits`	The class limits. Alternative to nclass.
`add`	logical. Whether a new graphic window is opened or the lines will be added into an exosting plot.
`se`	Logical. Use standard errors of means (se=TRUE, option (a) above) or class-specific standard deviations (se=FALSE, option (b) above).
`main`, `xlab`, `ylab`, `ylim`, `lwd`	Graphical parameters of the plot. ignored if `add=TRUE`.
`highlight`	The color for lines that do not cross the y-axis.

Details

The function first classifies the data in nclass classes of variable x so that each class has approximately equal number of observations. Then the class mean and deviation s is computed for each class, where s is either the standard error of the mean (if se=TRUE) or standard deviation (if se=FALSE). A vertical line is plotted at the middle of each class showing the class mean by a dot and lines of length 3.92*s. If the line does not cross the x- axis, then the highlight color is used in the line. With small number of observations (or lot of ties), the number of classes is decreased until each class includes the minimum of 2 observations.

Value

The function is usually used for its side effects (i.e., the plot). However, the values used in producing the plot are returned in a list of elements

x: the class middlepoint x values.
m: class-specific means of y.
s: class-specific standard deviations or standard errors of y (see details).
lb: lower ends of the class-specific lines.
ub: upper ends of the lines.

Author(s)

Lauri Mehtatalo <[email protected]>

References

Mehtatalo, Lauri and Lappi, Juha 2020. Biometry for Forestry and Environmental Data: with examples in R. New York: Chapman and Hall/CRC. 426 p. doi:10.1201/9780429173462

Examples

x<-seq(1,100,1)
y<-x+10*log(x)+rnorm(100,0,5)
fm1<-lm(y~x)
plot(x,resid(fm1))
mywhiskers(x,resid(fm1),se=FALSE,add=TRUE)
mywhiskers(x,resid(fm1),se=TRUE,lwd=2,add=TRUE)
abline(h=0)
x<-seq(1,100,1)
y<-x+10*log(x)+rnorm(100,0,5)
fm1<-lm(y~x)
plot(x,resid(fm1))
mywhiskers(x,resid(fm1),se=FALSE,add=TRUE)
mywhiskers(x,resid(fm1),se=TRUE,lwd=2,add=TRUE)
abline(h=0)

Solve a Nonlinear Equation Using Newton-Raphson algorithm.

Description

Solves an equation equation of form $f(x)=0$ , for scalar $x$ using the Newton-Raphson algorithm.

Usage

NR(init, fn, gr, crit = 6, range = c(-Inf, Inf))
NR(init, fn, gr, crit = 6, range = c(-Inf, Inf))

Arguments

`init`	Numeric scalar, The initial guess for $x$ .
`fn`	An R-function returning the scalar value of $f(x)$ , with $x$ as the only argument.
`gr`	An R-function returning the first derivative of $f(x)$ , with $x$ as the only argument.
`crit`	Convergence criteria. The upper limit for the absolute value of $f(x)$ at an accepted the solution.
`range`	A two-unit vector giving the upper and lower bounds for $x$ . The solution is searched from within this range.

Details

The function is a straightforward implementation of the well-known Newton-Raphson algorithm.

Value

A list of components

`par`	the value of $x$ in the solution
`crit`	the value of $f(x)$ at the solution

If estimation fails (no solution is found during 100000 iterations), both e lements of the solution are NA's.

Author(s)

Lauri Mehtatalo <[email protected]>

Examples

## Numerically solve Weibull shape for a stand
## where DGM=15cm, G=15m^2/ha and N=1000 trees per ha
func<-function(shape,G,N,DGM) {
##      print(G,DGM,N)
      val<-pi/(4*gamma(1-2/shape)*log(2)^(2/shape))-G/(N*DGM^2)
      val
      }

grad<-function(shape) {
      pi/4*(-1)*                                                            
      (gamma(1-2/shape)*log(2)^(2/shape))^(-2)*                             
      (gamma(1-2/shape)*digamma(1-2/shape)*2*shape^(-2)*log(2)^(2/shape)+   
      log(2)^(2/shape)*log(log(2))*(-2)*shape^(-2)*gamma(1-2/shape))         
      }
      
shape<-NR(5,fn=function(x) func(x,G=10000*15,1000,15),gr=grad,crit=10,range=c(2.1,Inf))$par
## Numerically solve Weibull shape for a stand
## where DGM=15cm, G=15m^2/ha and N=1000 trees per ha
func<-function(shape,G,N,DGM) {
##      print(G,DGM,N)
      val<-pi/(4*gamma(1-2/shape)*log(2)^(2/shape))-G/(N*DGM^2)
      val
      }

grad<-function(shape) {
      pi/4*(-1)*                                                            
      (gamma(1-2/shape)*log(2)^(2/shape))^(-2)*                             
      (gamma(1-2/shape)*digamma(1-2/shape)*2*shape^(-2)*log(2)^(2/shape)+   
      log(2)^(2/shape)*log(log(2))*(-2)*shape^(-2)*gamma(1-2/shape))         
      }
      
shape<-NR(5,fn=function(x) func(x,G=10000*15,1000,15),gr=grad,crit=10,range=c(2.1,Inf))$par

Solve a Systems of Nonlinear Equations Using the Newton's Method

Description

Solves a system of equations of form $f_1(x) = 0$ , $f_2(x) = 0$ ,..., $f_p(x) = 0$ for vector $x$ using the multidimensional version of the Newton-Raphson algorithm. The gradients are solved numerically within the function using R-function numericDeriv.

Usage

NRnum(init, fnlist, crit = 6, ...)
NRnum(init, fnlist, crit = 6, ...)

Arguments

`init`	vector of initial values for $x$ .
`fnlist`	list of R-functions for $f_1(x)$ , $f_2(x)$ , ..., $f_p(x)$ each function gets a vector-valued argument $x$ and returns a scalar value.
`crit`	Convergence criterion. Stop iteration when ( $\|f_1(x)\|+\|f_2(x)\|+...+\|f_p(x)\|<crit$ ).
`...`	Other arguments passed to the functions of `fnlist`

Value

A list of components

`par`	the value of vector $x$ in the solution
`crit`	the value of the convergence criterion at the solution

If estimation fails (no solution is found during 100 iterations), both elements of the solution are NA's.

Author(s)

Lauri Mehtatalo, <[email protected]>

Examples

# Moment-based recovery of Weibull parameters
mu<-14
mu2<-210
muf<-function(theta) theta[2]*gamma(1+1/theta[1])-mu
mu2f<-function(theta) theta[2]^2*gamma(1+2/theta[1])-mu2
functions<-list(muf,mu2f)
momrec<-NRnum(c(3,13),functions)
momrec$par
# Moment-based recovery of Weibull parameters
mu<-14
mu2<-210
muf<-function(theta) theta[2]*gamma(1+1/theta[1])-mu
mu2f<-function(theta) theta[2]^2*gamma(1+2/theta[1])-mu2
functions<-list(muf,mu2f)
momrec<-NRnum(c(3,13),functions)
momrec$par

Increment core data of Scots pine trees

Description

Growth ring measurements of 88 trees of a long-term thinning experiment on a naturally regenerated Scots pine stand in Eastern Finland.

Usage

data(patti)data(patti)

Format

A data frame with 3604 observations on the following 9 variables.

Plot: Sample plot id, a factor with 10 levels.
Tree: Tree id, a factor with 55 levels (same tree id may occur on different plots!).
SDClass: Thinning treatment, factor with 4 levels (1=Control, 2=Light, 3=Moderate, 4=Heavy).
Diam1986: Tree diameter in year 1986, just before the thinning.
Year: Calendar year of the ring.
CA: Current tree age in years.
RW: Ring width, $mm$
RD: Ring density, $g/cm^3$
RBA: Ring Basal area, $mm^2$

Details

Long-term thinning experiment on a naturally regenerated Scots pine stand in Eastern Finland. The experiment consists of 10 sample plots, in four different classes according to the post-thinning stand density. The plots were thinned in winter 1986-1987. In winter 2006–2007, 10 trees were felled from each plot. A radial 5mm by 5mm segment from pith to bark was cut from each tree at height 1.3 meter height. Ring widths from pith to bark were analyzed for each sample, using an ITRAX X-ray microdensitometer an post-processed to create ring widths from pith to bark were determined for each disc. The ring widths were further transformed to ring basal areas by assuming circular, growth rings. For 12 trees, ring widths could not be extracted. The data includes ring widths for a total of 88 trees between years 1991-2005.

References

Mehtatalo, Lauri and Lappi, Juha 2020. Biometry for Forestry and Environmental Data: with examples in R. New York: Chapman and Hall/CRC. 426 p. doi:10.1201/9780429173462.

Examples

data(afterthin)
par(mfcol=c(2,1),cex=0.7,mai=c(0.8,0.8,0.5,0.1))
linesplot(afterthin$CA,
		  afterthin$RBA,
		  group=afterthin$Plot:afterthin$Tree,
		  col.lin=as.numeric(afterthin$SDClass),cex=0,
		  xlab="Tree age",
		  ylab=expression("Ring basal area, "*mm^2))
		  
linesplot(afterthin$Year,
		  afterthin$RBA,
		  group=afterthin$Plot:afterthin$Tree,
		  col.lin=as.numeric(afterthin$SDClass),cex=0,
		  xlab="Year",
		  ylab=expression("Ring basal area, "*mm^2))
data(afterthin)
par(mfcol=c(2,1),cex=0.7,mai=c(0.8,0.8,0.5,0.1))
linesplot(afterthin$CA,
		  afterthin$RBA,
		  group=afterthin$Plot:afterthin$Tree,
		  col.lin=as.numeric(afterthin$SDClass),cex=0,
		  xlab="Tree age",
		  ylab=expression("Ring basal area, "*mm^2))
		  
linesplot(afterthin$Year,
		  afterthin$RBA,
		  group=afterthin$Plot:afterthin$Tree,
		  col.lin=as.numeric(afterthin$SDClass),cex=0,
		  xlab="Year",
		  ylab=expression("Ring basal area, "*mm^2))

The Percentile-based Distribution

Description

Density, distribution function, quantile function and random generation for the percentile-based distribution.

Usage

dPercbas(x, xi, F)
pPercbas(q, xi, F)
qPercbas(p, xi, F)
rPercbas(n, xi, F)
dPercbas(x, xi, F)
pPercbas(q, xi, F)
qPercbas(p, xi, F)
rPercbas(n, xi, F)

Arguments

`x`, `q`	vector of quantiles
`p`	vector of probabilitiies
`n`	number of observations. If `length(n) > 1`, the length is taken to be the number required.
`xi`	Strictly increasing vector of percentiles corresponding to the cumulative probabilities given in `F`. Of same length as `F`.
`F`	a `k`-length strictly increasing vector of cumulative probabilities, with `xi[1]=0` and `xi[k]=1`.

Details

The percentile-based distribution is defined by the quantiles xi that correspond to the cumulative probabilities given in F. The continuous distribution is obtained by linear interpolation of the cdf.

Value

dll gives the density, pll gives the distribution function, qll gives the quantile function, and rll generates random deviates.

The length of the result is determined by n for rPercbas, and by the length of x, q and p for the other functions.

Author(s)

Lauri Mehtatalo <[email protected]>

References

Borders B. E., Souter R. A., Bailey. R. L., and Ware, K. D. 1987. Percentile-based distributions characterize forest stand tables. Forest Science 33(2): 570-576.

Mehtatalo, L. 2005. Localizing a predicted diameter distribution using sample information. Forest Science 51(4): 292–302.

Mehtatalo, Lauri and Lappi, Juha 2020a. Biometry for Forestry and Environmental Data: with examples in R. New York: Chapman and Hall/CRC. 426 p. doi:10.1201/9780429173462

Mehtatalo, Lauri and Lappi, Juha 2020b. Biometry for Forestry and Environmental Data: with examples in R. Full Versions of The Web Examples. Available at http://www.biombook.org.

Examples

d0<-seq(0,30,0.01)
d<-c(5,7,10,11,11.7,13,15,19,22,24,25,28.5)
plot(d0,pPercbas(d0,d),type="l")
hist(rPercbas(1000,d),breaks=seq(0,30,1),freq=FALSE,ylim=c(0,0.15))
lines(d0,dPercbas(d0,d),col="red")
d0<-seq(0,30,0.01)
d<-c(5,7,10,11,11.7,13,15,19,22,24,25,28.5)
plot(d0,pPercbas(d0,d),type="l")
hist(rPercbas(1000,d),breaks=seq(0,30,1),freq=FALSE,ylim=c(0,0.15))
lines(d0,dPercbas(d0,d),col="red")

Sapling counts from sample plots of sapling stands in Finland.

Description

Grouped Norway Spruce regeneration establishment data.

Usage

data(plants)data(plants)

Format

A data frame with 1926 observations (fixed-area sample plots) from a total of 123 forest stands, with the following 8 variables.

spruces: The number of spruce saplings (both planted and natural)
stand: The stand id)
hdecid: The mean height of deciduous tree species
prepar: Site preparation method, categorical with 4 levels
stones: Binary indicator for stoniness
wet: Binary indicator for wetness

Details

The data are collected from 123 fixed-area sample plots with similar age of planted spruce saplings. The variables have been measured on fixed-area plots.

References

Miina, J. and Saksa, T. 2006. Predicting regeneration establishment in Norway spruce plantations using a multivariate multilevel model. New Forests 32: 265-283. doi:10.1007/s11056-006-9002-y

Mehtatalo, Lauri and Lappi, Juha 2020. Biometry for Forestry and Environmental Data: with examples in R. New York: Chapman and Hall/CRC. 426 p. doi:10.1201/9780429173462

Examples

data(plants)
library(lme4)
## Not run: 
glmm1<-glmer(spruces ~ (1|stand)+hdecid+as.factor(prepar)+as.factor(stones)+as.factor(wet), 
             family=poisson(), data=plants)

## End(Not run)
data(plants)
library(lme4)
## Not run: 
glmm1<-glmer(spruces ~ (1|stand)+hdecid+as.factor(prepar)+as.factor(stones)+as.factor(wet), 
             family=poisson(), data=plants)

## End(Not run)

Sapling counts from sample plots of sapling stands in Finland.

Description

Independent Norway Spruce regeneration establishment data.

Usage

data(plants2)data(plants2)

Format

A data frame with 123 observations on the following 8 variables.

planted: The number of planted spruce saplings on the plot
pines: The number of natural pine saplings
spruces: The number of natural spruce saplings
birches: The number natural birch saplings
othersp: The number of natural saplings of other species
hcrop: The mean height of crop species
hdecid: The mean height of deciduous tree species
sitetype: Site fertility class (small number indicates more fertile site)

Details

The data are collected from 123 fixed-area sample plots with similar age of planted spruce saplings. The number of saplings per species and the height of crop species (spruce and pine) and competing vegetation (birch and other broadleaved trees) has been recorded for all plots. The data includes one plot per forest stand.

References

Miina, J. and Saksa, T. 2006. Predicting regeneration establishment in Norway spruce plantations using a multivariate multilevel model. New Forests 32: 265-283. doi:10.1007/s11056-006-9002-y

Mehtatalo, Lauri and Lappi, Juha 2020. Biometry for Forestry and Environmental Data: with examples in R. New York: Chapman and Hall/CRC. 426 p. doi:10.1201/9780429173462

Examples

data(plants2)

glm1 <- glm(spruces ~ hdecid, family=quasipoisson(), data=plants2)
data(plants2)

glm1 <- glm(spruces ~ hdecid, family=quasipoisson(), data=plants2)

Diagnostic plot a Height-Diameter model residuals

Description

Plotting method for class hdmod

Usage

## S3 method for class 'hdmod'
plot(x, col.point = "blue", highlight = "red", standd = TRUE, 
           cex=1, corD=FALSE, ask=TRUE, ...)
## S3 method for class 'hdmod'
plot(x, col.point = "blue", highlight = "red", standd = TRUE, 
           cex=1, corD=FALSE, ask=TRUE, ...)

Arguments

`x`	A H-D model model fitted by `fithd`.
`col.point`	The color used for data points
`highlight`	The color used to higlight classes with mean sighnificantly different from zero
`standd`	Plot residuals against diameter standardized using the plot-specific mean and diameter (`standd=TRUE`) or against raw diameter (`standd=FALSE`)
`cex`	See `par`
`corD`	should predictions of random effects be plotted on mean diameter of the plot.
`ask`	ask before new plot.
`...`	Other arguments, currently ignored.

Details

The function makes residual plots on a fitted H-D model, which can be used to explore whether the fixed part satisfactorily models the shape of H-D models. The residuals are plotted on diameters standardized at plot level (dsd) or on raw diameters (d) according to argument standd. Here $dsd=(d-D)/SDD$ , where $d$ is tree diameter, $D$ and $SDD$ are the mean and standard deviation of diameters on the plot in question. Using plot-specific standardized diameter ensures that e.g., the medium-sized trees of the plot are always in the middle of the plot, which provides a better graph to explore the fit at the plot level in a dataset where the diameter range varies between plots.

Author(s)

Lauri Mehtatalo <[email protected]>

Examples

data(spati)

model<-fithd(spati$d,spati$h,spati$plot)

plot(model)
data(spati)

model<-fithd(spati$d,spati$h,spati$plot)

plot(model)

Individual tree volume functions for Finland

Description

Predict individual tree volumes using the functions of Laasasenaho(1982). The volume prediction can be based on tree diameter or tree diameter and height. The functions applying upper stem diameter have not (yet) been implemented.

Usage

predvol(species,d,h=0,model=1)
predvol(species,d,h=0,model=1)

Arguments

`species`	The vector of tree species. 1:Pine, 2:Spruce, 3: Silver birch. 4: Downy birch. Other codes than 1-4 are accepted but return NA as the volume prediction.
`d`	The vector of tree diameters at breast height (cm)
`h`	The vector of tree heights. Used only if model=2.
`model`	The model used. If model is 1, the prediction is based on tree diameter only. If model=2, then diameter and height are used.

Details

Vectors species, dbh and height should be either scalars or vectors of the same length so that each element corresponds to one individual tree.

Value

A vector of tree volumes (in liters).

Author(s)

Lauri Mehtatalo <[email protected]>

References

Laasasenaho, Jouko 1982. Taper curve and volume functions for pine, spruce and birch. Comm. Inst. For. Fenn 108: 1-74.

Examples

d<-c(15,18.3,29.3,22.4)
h<-c(13,18,22,19)
species<-c(1,1,1,3)
predvol(species,d,h,model=2)
predvol(species,d,model=1)
d<-c(15,18.3,29.3,22.4)
h<-c(13,18,22,19)
species<-c(1,1,1,3)
predvol(species,d,h,model=2)
predvol(species,d,model=1)

Normal QQ-plot of a fitted H-D model

Description

Produces a panel of graphs including the Normal qq-plot of a H-D model residuals and of the predicted random effects.

Usage

qqplotHD(model, startnew=TRUE)
qqplotHD(model, startnew=TRUE)

Arguments

`model`	A nonlinear H-D model model fitted by `fithd`.
`startnew`	Should a new plotting window be opened?

Details

The function extracts the residuals and the random effects of the fitted Height-Diameter model and produces a panel of plots including univariate Normal qq-plots of the model.

Author(s)

Lauri Mehtatalo <[email protected]>

Examples

data(spati)

model<-fithd(spati$d,spati$h,spati$plot)

qqplotHD(model)
data(spati)

model<-fithd(spati$d,spati$h,spati$plot)

qqplotHD(model)

Properties of sample quantiles from a tree population described by the percentile-based diameter distribution.

Description

Function qtree.moments finds the expected value and variance for $X_{r:n}$ ; the r:th smallest observation in an iid sample of size n from a population with a percentile-based distribution.

Function qtree.jointdens computes the bivariate pdf for two quantiles $(X_{r1:n},X_{r2:n})$ from the same sample, where $r1<r2$ .

Function qtree.exy approximates expected value of the product $X_{r1:n}X_{r2:n}$ , i.e. the integral of function $x_{r1:n}x_{r2:n}f_{r1:n,r2:n}(\bm x)$ over the two-dimensional range of $\bm x$ by computing for each percentile interval the function mean in a regular npts*npts grid and multiplying the mean by the area.

Function qtree.varcov returns the expected valuers, cumulative percentage values and the variance-covariance matrices that correspond to given sample quantiles and underlying percentile-based distribution of the population.

Function interpolate.D does a bilinear interpolation of the variance-covariance matrix of percentiles that correspond to values F of the cdf to values that correspond to values ppi.

Usage

qtree.moments(r,n,xi,F)
qtree.jointdens(x,r1,r2,n,xi,F)
qtree.exy(r1,r2,n,xi,F,npts=100) 
qtree.varcov(obs,xi,F)
interpolate.D(D,ppi,F) 
qtree.moments(r,n,xi,F)
qtree.jointdens(x,r1,r2,n,xi,F)
qtree.exy(r1,r2,n,xi,F,npts=100) 
qtree.varcov(obs,xi,F)
interpolate.D(D,ppi,F)

Arguments

`r`, `r1`, `r2`	The ranks of the sample order statistics. `r=1` means the smallest, `r=n` the largest.
`n`	The sample size
`xi`	The percentiles that specify the cdf in increasing order. The first element should be the population minimum and the last element should be the population maximum. A vector of same length as `F`
`F`	The values of the cdf that correspond to the percentiles of `xi`. The first elements should be 0 and the last 1.
`x`	a matrix with two columns that gives the x-values for which the joint density is computed in `qtree.jointdens`.
`npts`	The number of regularly placed points that is used in the integral approximation of $E(X_{r1:n}X_{r2:n})$ for each percentile interval in function `exy`.
`obs`	A data frame of observed sample quantiles, possibly from several plots. The data frame should include (at least) columns `r` (the ranks), `n` (sample size), `plot` (plot id) and `d` (observed diameter). The rows should be ordered by `r` within each plot, and all observations from same plot should follow each other.
`D`	The variance-covariance matrix of the residual errors (plot effects) of percentile models. The number of columns and rows should equal to the length of `F` and `xi`.
`ppi`	The values of cdf for which the covariances needs to be interpolated in `interpolate.D`.

Value

Function qtree.moments returns a list with elements

`mu`	The expected value of $X_{r:n}$ .
`sigma2`	The variance of $X_{r:n}$ .
`x`, `y`	y gives the values of the pdf of $X_{r:n}$ for values given in x for plotting purposes. Try `plot(sol$x,sol$y,type="l")`.

Function qtree.jointdens returns a vector with length equal to the nrow(x), including the values of the joint pdf of $({X_{r1:n}},X_{r2:n})$ in these points.

Function qtree.exy returns a scalar, the approximate of $E(X_{r1:n}X_{r2:n})$ .

Function qtree.varcov returns a list with elements

`obs`	The original input data frame, augmented with the expected values in column `Ed` and the corresponding values of the cdf of $X$ in column `pEd`.
`R`	The variance-covariance matrix of the sample quantiles.

Function interpolate.D returns a list with elements

`D`	The original variance-covariance matrix, augmented with the variances and covariances that correspond to the cdf values `ppi`.
`F`	The values of cdf that correspond to the augmented matrix `D`.
`D1`	The variance-covariance matrix of the percentiles that correspond to the cdf values given in `ppi`
`D2`	The covariance matrix between the percentiles that correspond to `ppi` and `F`

Author(s)

Lauri Mehtatalo <[email protected]>

References

Mehtatalo, L. 2005. Localizing a predicted diameter distribution using sample information. Forest Science 51(4): 292–302.

Mehtatalo, Lauri and Lappi, Juha 2020a. Biometry for Forestry and Environmental Data: with examples in R. New York: Chapman and Hall/CRC. 426 p. doi:10.1201/9780429173462

Mehtatalo, Lauri and Lappi, Juha 2020b. Biometry for Forestry and Environmental Data: with examples in R. Full Versions of The Web Examples. Available at http://www.biombook.org.

Examples

F<-c(0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,0.95,1)

# Predictions of logarithmic percentiles
xi<-c(1.638,2.352,2.646,2.792,2.91,2.996,3.079,3.151,3.234,3.349,3.417,3.593)

# The variance of their prediction errors
D<-matrix(c(0.161652909,0.050118692,0.022268974,0.010707222,0.006888751,0,
   0.000209963,-0.002739361,-0.005478838,-0.00655718,-0.006718843,-0.009819052,
   0.050118692,0.074627668,0.03492943,0.01564454,0.008771398,0,
   -0.002691651,-0.005102312,-0.007290366,-0.008136685,-0.00817717,-0.009026883,
   0.022268974,0.03492943,0.029281808,0.014958206,0.009351904,0,
   -0.002646641,-0.003949305,-0.00592412,-0.006556639,-0.006993025,-0.007742731,
   0.010707222,0.01564454,0.014958206,0.014182608,0.009328299,0,
   -0.001525745,-0.002448765,-0.003571811,-0.004470387,-0.004791053,-0.005410252,
   0.006888751,0.008771398,0.009351904,0.009328299,0.009799233,0,
   -0.000925308,-0.001331631,-0.002491679,-0.003277911,-0.003514961,-0.003663479,
   rep(0,12),
   0.000209963,-0.002691651,-0.002646641,-0.001525745,-0.000925308,0,
   0.003186033,0.003014887,0.002961818,0.003112953,0.003050486,0.002810937,
   -0.002739361,-0.005102312,-0.003949305,-0.002448765,-0.001331631,0,
   0.003014887,0.00592428,0.005843888,0.005793879,0.005971638,0.006247869,
   -0.005478838,-0.007290366,-0.00592412,-0.003571811,-0.002491679,0,
   0.002961818,0.005843888,0.00868157,0.008348973,0.008368812,0.008633202,
   -0.00655718,-0.008136685,-0.006556639,-0.004470387,-0.003277911,0,
   0.003112953,0.005793879,0.008348973,0.011040791,0.010962609,0.010906917,
   -0.006718843,-0.00817717,-0.006993025,-0.004791053,-0.003514961,0,
   0.003050486,0.005971638,0.008368812,0.010962609,0.013546621,0.013753718,
   -0.009819052,-0.009026883,-0.007742731,-0.005410252,-0.003663479,0,
   0.002810937,0.006247869,0.008633202,0.010906917,0.013753718,0.02496596),ncol=12)

# observed tree data, 5 trees from 2 plots
obs<-data.frame(r=c(1,3,6,1,2),n=c(7,7,7,9,9),plot=c(1,1,1,2,2),d=c(10,11,27,8,12))

# See Example 11.33 in Mehtatalo and Lappi 2020b
qtrees<-qtree.varcov(obs,xi,F)
obs<-qtrees$obs
mustar<-obs$Ed
ystar<-log(obs$d)
R<-qtrees$R
Dtayd<-interpolate.D(D,obs$pEd)
F<-c(0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,0.95,1)

# Predictions of logarithmic percentiles
xi<-c(1.638,2.352,2.646,2.792,2.91,2.996,3.079,3.151,3.234,3.349,3.417,3.593)

# The variance of their prediction errors
D<-matrix(c(0.161652909,0.050118692,0.022268974,0.010707222,0.006888751,0,
   0.000209963,-0.002739361,-0.005478838,-0.00655718,-0.006718843,-0.009819052,
   0.050118692,0.074627668,0.03492943,0.01564454,0.008771398,0,
   -0.002691651,-0.005102312,-0.007290366,-0.008136685,-0.00817717,-0.009026883,
   0.022268974,0.03492943,0.029281808,0.014958206,0.009351904,0,
   -0.002646641,-0.003949305,-0.00592412,-0.006556639,-0.006993025,-0.007742731,
   0.010707222,0.01564454,0.014958206,0.014182608,0.009328299,0,
   -0.001525745,-0.002448765,-0.003571811,-0.004470387,-0.004791053,-0.005410252,
   0.006888751,0.008771398,0.009351904,0.009328299,0.009799233,0,
   -0.000925308,-0.001331631,-0.002491679,-0.003277911,-0.003514961,-0.003663479,
   rep(0,12),
   0.000209963,-0.002691651,-0.002646641,-0.001525745,-0.000925308,0,
   0.003186033,0.003014887,0.002961818,0.003112953,0.003050486,0.002810937,
   -0.002739361,-0.005102312,-0.003949305,-0.002448765,-0.001331631,0,
   0.003014887,0.00592428,0.005843888,0.005793879,0.005971638,0.006247869,
   -0.005478838,-0.007290366,-0.00592412,-0.003571811,-0.002491679,0,
   0.002961818,0.005843888,0.00868157,0.008348973,0.008368812,0.008633202,
   -0.00655718,-0.008136685,-0.006556639,-0.004470387,-0.003277911,0,
   0.003112953,0.005793879,0.008348973,0.011040791,0.010962609,0.010906917,
   -0.006718843,-0.00817717,-0.006993025,-0.004791053,-0.003514961,0,
   0.003050486,0.005971638,0.008368812,0.010962609,0.013546621,0.013753718,
   -0.009819052,-0.009026883,-0.007742731,-0.005410252,-0.003663479,0,
   0.002810937,0.006247869,0.008633202,0.010906917,0.013753718,0.02496596),ncol=12)

# observed tree data, 5 trees from 2 plots
obs<-data.frame(r=c(1,3,6,1,2),n=c(7,7,7,9,9),plot=c(1,1,1,2,2),d=c(10,11,27,8,12))

# See Example 11.33 in Mehtatalo and Lappi 2020b
qtrees<-qtree.varcov(obs,xi,F)
obs<-qtrees$obs
mustar<-obs$Ed
ystar<-log(obs$d)
R<-qtrees$R
Dtayd<-interpolate.D(D,obs$pEd)

Recovery of Weibull parameters of tree diameter distribution using measured stand characteristics

Description

The function finds such parameters shape and scale of the Weibull diameter distribution that yield the given basal area, number of stems and weighted/unweighted mean/median diameter. Weibull function can be assumed either as the unweighted or basal-area weighted distribution.

Usage

recweib(G, N, D, Dtype, init=NA, trace=FALSE, weight=0,minshape=0.01)
func.recweib1(lshape, G, N, D, Dtype, trace=FALSE,minshape=0.01)
func.recweib2(lshape, G, N, D, Dtype, trace=FALSE,minshape=0.01)
recweib(G, N, D, Dtype, init=NA, trace=FALSE, weight=0,minshape=0.01)
func.recweib1(lshape, G, N, D, Dtype, trace=FALSE,minshape=0.01)
func.recweib2(lshape, G, N, D, Dtype, trace=FALSE,minshape=0.01)

Arguments

`G`	The basal area in $m^2/ha$ , scalar.
`N`	The number of stems per ha, scalar.
`D`	Either A: The arithmetic mean diameter, B: The basal-area weighted mean diameter, C: median diameter or D: The basal-area weighted median diameter of the stand, cm.
`Dtype`	One of characters "A", "B", "C", "D", indicating which type of mean diameter was given in argument D.
`init`	The initial guess for the shape parameter (scalar). If not given, a simple model (see Siipilehto and Mehtatalo 2013, appendix) is used to compute the initial guess for the unweighted case; value 4 is used as default in the basal-area weighted case.
`trace`	if TRUE, some output on the convergence of the algorithm is printed on the screen.
`weight`	if weight=0 (the default), Weibull function is assumed as the unweighted density. If weight=2, weibull function is assumed as the basal-area weighted density. Weight is also the theoretical infimum of the shape parameter.
`lshape`	logarithmic shape parameter, (log(shape+minshape)).
`minshape`	Minimum difference of the shape parameter and its theoretical infimum.

Details

The recovery is based on the solution of the equation DQMW^2(shape,scale(D,shape))-DQM^2= 0, where DQMW(shape, scale(D,shape)) expresses the DQM of the assumed Weibull distribution for the given value of the shape parameter and using the scale parameter that corresponds to the given combination of the shape parameter and the mean/median diameter given in D. The function which is set to zero is implemented in functions func.recweib1 (unweighted case) and func.recweib2 (ba-weighted case). The Gauss-Newton method implemented in NRnum is used for solving the equation.

Value

A list of components

`shape`, `scale`	The value of the shape and scale parameters at the solution.
`G`, `N`, `D`, `Dtype`	The input arguments.
`val`	The value of the equation DQMW^2(shape,scale(D,shape))-DQM^2 at the solution

Author(s)

Lauri Mehtatalo <[email protected]> and Jouni Siipilehto

References

Siipilehto, J. and Mehtatalo, L. 2013. Parameter recovery vs. parameter prediction for the Weibull distribution validated for Scots pine stands in Finland. Silva Fennica 47(4), article id 1057. doi:10.14214/sf.1057

Mehtatalo, Lauri and Lappi, Juha 2020a. Biometry for Forestry and Environmental Data: with examples in R. New York: Chapman and Hall/CRC. 426 p. doi:10.1201/9780429173462

Mehtatalo, Lauri and Lappi, Juha 2020b. Biometry for Forestry and Environmental Data: with examples in R. Full Versions of The Web Examples. Available at http://www.biombook.org.

Examples

# Demonstration with 3 example stands.
# Example stand 1. Uneven-aged stand in Finland (Vesijako, Kailankulma, stand no 1):
G<-17.0
N<-1844
D<-7.9
DG<-19.6
DM<-8.1
DGM<-19.1
recweib(G,N,D,"A")            #  1.066123 8.099707
recweib(G,N,DG,"B")           # 1.19316  8.799652
recweib(G,N,DM,"C")           # 1.601795 10.18257
recweib(G,N,DGM,"D")          # 1.095979 8.280063
recweib(G,N,D,"A",weight=2)   # 2.590354 21.66398
recweib(G,N,DG,"B",weight=2)  # 2.563892 22.07575
recweib(G,N,DM,"C",weight=2)  # 2.998385 17.74291
recweib(G,N,DGM,"D",weight=2) # 2.566621 22.03183



# Example 2. Even aged stand in Finland (see Siipilehto & Mehtatalo, Fig 2):
G_ha<-9.6
N_ha<-949
D<-11.0
DG<-12.3
DM<-11.1
DGM<-12.4  
recweib(G_ha,N_ha,D,"A") # 4.465673 12.05919
recweib(G_ha,N_ha,DG,"B") # 4.463991 12.05912
recweib(G_ha,N_ha,DM,"C")  # 4.410773 12.05949
recweib(G_ha,N_ha,DGM,"D") # 4.448272 12.05924


# Example 3. Assumed peaked even aged stand (see Siipilehto & Mehtatalo, Fig 1):
G_ha<-10.0
N_ha<-1300
D<-9.89
DG<-10.0
DM<-9.89
DGM<-10.0  
recweib(G_ha,N_ha,D,"A")  #  34.542 10.04978
recweib(G_ha,N_ha,DG,"B") # 14.23261 10.22781
recweib(G_ha,N_ha,DM,"C") # 6.708882 10.44448
recweib(G_ha,N_ha,DGM,"D") # 24.45228 10.10607
# Demonstration with 3 example stands.
# Example stand 1. Uneven-aged stand in Finland (Vesijako, Kailankulma, stand no 1):
G<-17.0
N<-1844
D<-7.9
DG<-19.6
DM<-8.1
DGM<-19.1
recweib(G,N,D,"A")            #  1.066123 8.099707
recweib(G,N,DG,"B")           # 1.19316  8.799652
recweib(G,N,DM,"C")           # 1.601795 10.18257
recweib(G,N,DGM,"D")          # 1.095979 8.280063
recweib(G,N,D,"A",weight=2)   # 2.590354 21.66398
recweib(G,N,DG,"B",weight=2)  # 2.563892 22.07575
recweib(G,N,DM,"C",weight=2)  # 2.998385 17.74291
recweib(G,N,DGM,"D",weight=2) # 2.566621 22.03183



# Example 2. Even aged stand in Finland (see Siipilehto & Mehtatalo, Fig 2):
G_ha<-9.6
N_ha<-949
D<-11.0
DG<-12.3
DM<-11.1
DGM<-12.4  
recweib(G_ha,N_ha,D,"A") # 4.465673 12.05919
recweib(G_ha,N_ha,DG,"B") # 4.463991 12.05912
recweib(G_ha,N_ha,DM,"C")  # 4.410773 12.05949
recweib(G_ha,N_ha,DGM,"D") # 4.448272 12.05924


# Example 3. Assumed peaked even aged stand (see Siipilehto & Mehtatalo, Fig 1):
G_ha<-10.0
N_ha<-1300
D<-9.89
DG<-10.0
DM<-9.89
DGM<-10.0  
recweib(G_ha,N_ha,D,"A")  #  34.542 10.04978
recweib(G_ha,N_ha,DG,"B") # 14.23261 10.22781
recweib(G_ha,N_ha,DM,"C") # 6.708882 10.44448
recweib(G_ha,N_ha,DGM,"D") # 24.45228 10.10607

The Weibull scale parameter for the given mean/median diameter and shape parameter.

Description

The function finds such scale parameter of the Weibull distribution that yields the given mean/median diameter. Function scaleDMean is used for arithmetic mean, scaleDGMean for the mean of basal-area weighted distribution, scaleDMed for median and scaleDGMed for the median of the basal-area weighted diameter distribution. Functions with number 1 in the name use Weibull functions as the unweighted density and functions with value 2 in the name use Weibull function as the basal-area weighted density.

The functions are used in the recovery of Weibull parameters using function recweib.

Usage

scaleDMean1(D,shape)
scaleDGMean1(D,shape)
scaleDMed1(D,shape)
scaleDGMed1(D,shape)
scaleDMean2(D,shape)
scaleDGMean2(D,shape)
scaleDMed2(D,shape)
scaleDGMed2(D,shape)
scaleDMean1(D,shape)
scaleDGMean1(D,shape)
scaleDMed1(D,shape)
scaleDGMed1(D,shape)
scaleDMean2(D,shape)
scaleDGMean2(D,shape)
scaleDMed2(D,shape)
scaleDGMed2(D,shape)

Arguments

`D`	The diameter
`shape`	The Weibull shape parameter

Value

scale

The value of the Weibull scale parameter.

Author(s)

Lauri Mehtatalo <[email protected]> and Jouni Siipilehto

References

Examples

scaleDMean1(15,3)
scaleDGMean1(15,3)
scaleDMed2(15,3)
scaleDGMed2(15,3)
scaleDMean1(15,3)
scaleDGMean1(15,3)
scaleDMed2(15,3)
scaleDGMed2(15,3)

Raw sample plot data of Scots pine in Ilomantsi, Finland.

Description

A dataset of Scots pine growth. The trees were collected on 56 fixed-area sample plots. The data includes no remeasurements. The growth data are based on measurements of increment borer chips.

Usage

data(spati)data(spati)

Format

A data frame with 9913 observations on the following variables.

plot: A unique sample plot id.
X: Plot size in X-direction, meters
Y: Plot size in y-direction, meters
N: Stand density, trees per ha
G: Basal area, $m^2/ha$
V: Plot volume, $m^3/ha$
Dg: Basal-area weighted mean diameter, cm
Hg: Height of basal area median diameter tree, m
Tg: Age of basal area median tree, yr
Hdom: Dominant height, m
maos: percentage of Scots pines of the total volume
kuos: percentage of Norway spruces of the total volume
kanro: A unique sample plot id (same as plot).
puunro: Tree id within plot.
pl: tree species. 1=Scots Pine
xk: x- coordinates of trees within plot
yk: y- coordinates of trees within plot
d: Tree diameter at breast height (1.3 meters above the ground) in cm.
h: Tree height, m.
t: Tree age, years
dk: Tree diameter at stump height, cm. there seems to be some unclear issues.
X2b: Double bark thickness, mm
id1: Tree diameter growth within the 5 year period prior to the measurement. Missing data coded as -1.
id2: Tree diameter growth within the period 6-10 years prior to the measurement. Missing data coded as -1.

Author(s)

The data were collected by Timo Pukkala.

References

Pukkala, T. 1989. Prediction of tree diameter and height in a Scots pine stand as a function of the spatial pattern of trees. Silva Fennica 23(2): 83-99. doi:10.14214/sf.a15532

Mehtatalo, Lauri and Lappi, Juha 2020. Biometry for Forestry and Environmental Data: with examples in R. New York: Chapman and Hall/CRC. 426 p. doi:10.1201/9780429173462

Heights and diameters of Scots pine trees in Ilomantsi, Finland.

Description

A dataset of Scots pine tree heights and diameters. The trees were collected on 56 fixed-area sample plots. This is a subset of the larger data set spati.

Usage

data(spati2)data(spati2)

Format

A data frame with 1678 observations on the following 3 variables.

plot: A unique sample plot id.
d: Tree diameter at breast height (1.3 meters above the ground) in cm.
h: Tree height, m.
n: The total number of trees on the plot.
dvar: The variance of tree diameters on the plot.
dmean: The mean of tree diameters on the plot.

Author(s)

The data were collected by Timo Pukkala.

References

Pukkala, T. 1989. Prediction of tree diameter and height in a Scots pine stand as a function of the spatial pattern of trees. Silva Fennica 23(2): 83-99. doi:10.14214/sf.a15532

Mehtatalo, Lauri and Lappi, Juha 2020. Biometry for Forestry and Environmental Data: with examples in R. New York: Chapman and Hall/CRC. 426 p. doi:10.1201/9780429173462

Examples

data(spati2)
fithd(spati2$d,spati2$h,spati2$plot)
data(spati2)
fithd(spati2$d,spati2$h,spati2$plot)

Productivity of stump lifting machines.

Description

The productivity of stump lifting machines on three Norway Spruce (Picea Abies) clearcut areas (sites). Stumps are lifted for use as bioenergy. The data were collected from three sites in Central Finland.

Usage

data(stumplift)data(stumplift)

Format

A data frame with 485 observations on the following 5 variables.

Stump: A unique stump id based on the order of processing. The successive numbers are usually close to each other in the clearcut area, but nearby trees do not necessarily have small difference in stump id.
Machine: The machine/clearcut/dirver combination. A factor with three levels.
Diameter: Stump diameter, cm.
Time: Processing time, seconds.
Productivity: Productivity, $m^3$ /effective working hour

Details

Each site was operated with different machine and driver so that the effect of site, machine and driver cannot be separated. The volume of each stump was estimated using the function of Laitila (2008), based on the stump diameter. A work system study was conducted to measure the processing time (seconds) and productivity ( $m^3$ /hour) for each stump.

References

Teijo Palander, Kalle Karha, Lauri Mehtatalo 2016. Applying polynomial regression modeling to productivity analysis of sustainable stump harvesting. Scandinavian Journal of Forest Reseach. doi:10.1080/02827581.2016.1238957

Teijo Palander, Janne Smolander, Kalle Karha, 2015. Work system study of three stump-lifting devices in Finland. Scandinavian Journal of Forest Research 30(6) 558-567, doi:10.1080/02827581.2015.1027731

Mehtatalo, Lauri and Lappi, Juha 2020. Biometry for Forestry and Environmental Data: with examples in R. New York: Chapman and Hall/CRC. 426 p. doi:10.1201/9780429173462

Examples

data(stumplift)
library(nlme)

modConstPow<-gls(Productivity~Machine+Machine*I((Diameter-70)^2),
                 data=stumplift,
                 weights=varPower(),
                 corr=corAR1(form=~Stump|Machine))
           
data(stumplift)
library(nlme)

modConstPow<-gls(Productivity~Machine+Machine*I((Diameter-70)^2),
                 data=stumplift,
                 weights=varPower(),
                 corr=corAR1(form=~Stump|Machine))

Effect of thinning on individual tree growth for 62 trees.

Description

Time series of estimated effect of thinning on the annual basal area growth of 62 Scots pine trees from three different thinning intensities.

Usage

data(thefdata)data(thefdata)

Format

A data frame with 1238 observations on the following 7 variables.

Plot: Sample plot id, a factor with 10 levels.
Tree: Tree id, a factor with 55 levels (same tree id may occur on different plots!).
Year: Calendar year of the ring.
SDClass: Thinning treatment, factor with 4 levels (1=Control, 2=Light, 3=Moderate, 4=Heavy).
CA: Current tree age in years.
Diam1986: Tree diameter in 1986 (just before the thinning).
ThEf: The thinning effect in annual basal area growth, mm^2.

Details

The data are based on a long-term thinning experiment on a naturally regenerated Scots pine stand in Eastern Finland, see documentation of afterthin and patti for details. The experiment consists of 10 sample plots, in four different classes (including an unthinned control) according to the post-thinning stand density. The plots were thinned in winter 1986-1987. The undisturbed growth for each tree was predicted using a linear mixed-effect model with crossed year and tree effect, which was fitted to a data set including the control treatment for all calendar years and other trees until the year of thinning. The thinning effect was computed as the difference of observed growth and the undisturbed predicted growth. For details about the procedure used in extracting the thinning effect, see Example 6.6 in Mehtatalo and Lappi 2020b and for nonlinear modeling of this data, see Chapter 7 of Mehtatalo and Lappi 2020a. The data includes only the observations after the thinning year for the tree thinned treatments (i.e., control is excluded).

References

Mehtatalo, Lauri and Lappi, Juha 2020a. Biometry for Forestry and Environmental Data: with examples in R. New York: Chapman and Hall/CRC. 426 p. doi:10.1201/9780429173462

Mehtatalo, Lauri and Lappi, Juha 2020b. Biometry for Forestry and Environmental Data: with examples in R. Full Versions of The Web Examples. Available at http://www.biombook.org.

Examples

data(thefdata)
linesplot(thefdata$Year,thefdata$ThEf, thefdata$Tree,col.lin=thefdata$SDClass)
data(thefdata)
linesplot(thefdata$Year,thefdata$ThEf, thefdata$Tree,col.lin=thefdata$SDClass)

Effect of thinning on individual tree growth

Description

A time series of estimated effect of thinning on the annual basal area growth of a Scots pine tree.

Usage

data(thinning)data(thinning)

Format

A data frame with 23 observations on the following 3 variables.

TreeID: Tree ID, 3_3 for all observations in this data.
Year: Calendar year (1983-2005).
ThEff: Estimated effect of thinning on the annual basal area growth in mm^2.

Details

The thinning took place between years 1986 and 1987. For details about the original measurements, see the documentation of data set patti, afterthin. For details about the procedure used in extracting the thinning effect, see Example 6.6 in Mehtatalo and Lappi 2020b and for nonlinear modeling of this data, see Chapter 7 of Mehtatalo and Lappi 2020a.

References

Mehtatalo, Lauri and Lappi, Juha 2020a. Biometry for Forestry and Environmental Data: with examples in R. New York: Chapman and Hall/CRC. 426 p. doi:10.1201/9780429173462

Mehtatalo, Lauri and Lappi, Juha 2020b. Biometry for Forestry and Environmental Data: with examples in R. Full Versions of The Web Examples. Available at http://www.biombook.org.

Examples

data(thinning)
plot(thinning$Year,thinning$ThEff,type="l")
data(thinning)
plot(thinning$Year,thinning$ThEff,type="l")

Individual tree volume modeling data for 8508 pine, spruce and birch trees in Finland

Description

Measurements of tree volumes from 8508 individual Scots pine (4066 trees), Norway spruce (2966 trees) and birch (1476 trees) trees in Finland. A total of 5053 trees were measured in 1970's by climbing to the standing trees, 1534 trees by falling the trees in 1990's, and 1921 trees by using terrestrial laser scanning in 2010's. The trees have been selected from among the sample trees of national forest inventories. In addition to tree volumes, diameters, heights and upper diameters, the data set includes variables describing the grouping of the trees to stands and plots and information about the site characteristics.

Usage

data(treevol)data(treevol)

Format

A data frame with 8508 observations on the following 11 variables.

stand: Stand id (numeric)
plot: Plot id (nested within stand)
v: Tree volume, dm^3
dbh: Tree diameter at 1.3m height above the stump height, cm
h: Tree height, m
d6: Tree diameter at 6m height above the stump height, cm
dataset: Categorical indicator for the dataset, (see the general description).
temp_sum: Numeric temperature sum of the site. A high value indicates warm site.
species: Categorical with three levels, tree species.
soil: Categorical with two levels: mineral soil or peatland.
region: Categorical with four levels: region.

Details

The data combines three datasets: 1) The climbed data including 2326 Scots pine (Pinus sylvestris L.), 1864 Norway spruce (Picea abies L. Karst.) and 863 silver (Betula pendula Roth) and downy (Betula pubescens Ehrh.) birches collected 1968-1972, and 2) the felled data with 797 pines, 479 spruces and 258 birches collected 1988-2001, and the scanned data with 943 Scots pine, 623 Norway spruce and 355 downy or silver birch collected 2017-2018. For the two older dat sets, diameter measurements are based on two perdendicular measurements of the stem at regularly spaced height. In the scanned data, the diemater measurements are based on TLS but heights are measured manually. For more details concerning the climbed data, see Laasasenaho (1982), and for the scanned data see Pitkanen et al. (2021) and for the felled data see (Korhonen & Maltamo 1990). In the climbed and felled data, the volumes of the trees were predicted using a natural cubic spline (splinefun, R Core Team 2021) to interpolate the taper curve between measured diameters at given relative heights. The volumes were obtained as a slolid of revolution based on the interpolated taper curve. In scanned data, a smoothing spline (smooth.spline, R core team 2021) function was used instead of an interpolating spline. The pooled data has been used in Kangas et al. (XXXX).

References

Kangas A., Pitkanen T., Mehtatalo L., Heikkinen J. (xxxx) Mixed linear and non-linear tree volume models with regionally varying parameters

Korhonen, K. and Maltamo, M. 1990. Mannyn maanpaallisten osien kuivamassat Etela-Suomessa. Metsantutkimuslaitoksen tiedonantoja 371.

Laasasenaho, J. (1982). Taper curve and volume functions for pine, spruce and birch. Communicationes Instituti Forestalia Fennica 108: 1-74.

Pitkanen, T.P., Raumonen, P., Liang, X., Lehtomaki, M. and Kangas, A.S. 2020. Improving TLS-based stem volume estimates by field measurements. Computers and Electronics in Agriculture and Forestry 180, 105882.

Examples

## Not run: 
data(treevol)
        
treevol$formfactor<-treevol$v/volvff(treevol$dbh,treevol$h,logita=100,lambda=log(0.2))
treevol$logitff<-log((treevol$formfactor)/(1-(treevol$formfactor)))

ptrees<-treevol[treevol$species=="pine",]
mod.init<-lm(logitff~I(1/h)+h+dbh+I(h*dbh)+I(1/(h*dbh))+
                     dataset+dataset:dbh+dataset:h,data=ptrees)

mod<-nlme(v~volvff(dbh,h,logita=logita,lambda=lambda),
         fixed=list(logita~I(1/h)+h+dbh+I(h*dbh)+I(1/(h*dbh))+dataset+
                           dataset:dbh+dataset:h+soil+temp_sum,
                           lambda~1),
         random=logita~1|stand/plot,
         start=c(coef(mod.init),rep(0,2),log(0.2)),
         data=ptrees,
         weights=varComb(varIdent(form=~1 |dataset),varPower()),
         method="ML",
         control=list(msVerbose=TRUE),
         verbose=TRUE)

## End(Not run)                    
## Not run: 
data(treevol)
        
treevol$formfactor<-treevol$v/volvff(treevol$dbh,treevol$h,logita=100,lambda=log(0.2))
treevol$logitff<-log((treevol$formfactor)/(1-(treevol$formfactor)))

ptrees<-treevol[treevol$species=="pine",]
mod.init<-lm(logitff~I(1/h)+h+dbh+I(h*dbh)+I(1/(h*dbh))+
                     dataset+dataset:dbh+dataset:h,data=ptrees)

mod<-nlme(v~volvff(dbh,h,logita=logita,lambda=lambda),
         fixed=list(logita~I(1/h)+h+dbh+I(h*dbh)+I(1/(h*dbh))+dataset+
                           dataset:dbh+dataset:h+soil+temp_sum,
                           lambda~1),
         random=logita~1|stand/plot,
         start=c(coef(mod.init),rep(0,2),log(0.2)),
         data=ptrees,
         weights=varComb(varIdent(form=~1 |dataset),varPower()),
         method="ML",
         control=list(msVerbose=TRUE),
         verbose=TRUE)

## End(Not run)

Solve a simple equation using a step halving algorithm.

Description

Solves equations of form $f(x)=0$ , for scalar $x$ $(l<=x<=u)$ using a simple step halving algorithm, where $f(x)$ is a monotonic continuous function. Initial finite upper and lower bounds for x are required. The algorithm first computes $f$ for $x=u$ and $x=l$ . If the sign was different then another call is performed at the midpoint $x=(u+l)/2$ , and the midpoint is taken as a new upper or lower bound, according to the location of sign change. The upper or lower bound are repeatedly updated until the absolute value of f at the midpoint is below a specified criteria.

Usage

updown(l, u, fn, crit = 6)
updown(l, u, fn, crit = 6)

Arguments

`l`	The initial lower bound
`u`	The initial upper bound
`fn`	R-function for $f(x)$
`crit`	The convergence criteria (Maximum accepted value of f at the solution is $10^{-\code{crit}}$ ).

Value

A scalar giving the value of $x$ at the solution. If the sign did not change between l and u, NA is returned.

Warning

May lead to infinite loop for non-continuous functions. Works only with monotonic functions.

Author(s)

Lauri Mehtatalo <[email protected]>

Examples

## Compute the median of Weibull distibution
fn<-function(x) pweibull(x,5,15)-0.5
updown(1,50,fn)
## Compute the median of Weibull distibution
fn<-function(x) pweibull(x,5,15)-0.5
updown(1,50,fn)

The Variable Form-Factor Volume Model

Description

An R-function for the variable form-factor volume model and a function for computing bias-corrected volumes augmented with parameter uncertainty from a model fitted using nlme.

Usage

volvff(dbh,h,theta=NA,logita=NA,lambda=NA)
predvff(data,mod,p=0.05,varMethod="taylor",biasCorr="none",nrep=500)
volvff(dbh,h,theta=NA,logita=NA,lambda=NA)
predvff(data,mod,p=0.05,varMethod="taylor",biasCorr="none",nrep=500)

Arguments

`dbh`	vector of individual tree diameters, cm
`h`	vector of individual tree heights (m), of same length as `dbh`.
`theta`, `logita`, `lambda`	The parameters `logita` and `lambda` can be given either in a two-column matrix where the number of rows equals to the length of `dbh`, or separately in objects `logita` and `lambda`, which can be scalars or vectors of same length as `dbh`. If `theta` is not `NA`, argumets `logita` and `lambda` are ignored.
`data`	A data set including variables `dbh` and `h` and all variables that are used as second-order predictors in `mod`.
`mod`	A variable form-factor model fitted using `nlme`. Any second-order predictors can be used for parameters `logita` and `lambda`. The bias-correction works properly only if the model includes random intercepts in `logita`.
`p`	The probability used in constructiong the confidence intervals for tree-level volumes and total volume. Symmetric $100(1-p)$ % confidence intervals are returned.
`varMethod`	Either `"taylor"`(the default) or `"simul"`. Specifies whether parameter uncertainty is approximated by using a linearized model based on the first-order Taylor approximation or by using Monte Carlo simulation.
`biasCorr`	Either `"none"`(the default), `"integrate"` or `"twopoint"`. Specifies whether the bias correction is not done (`"none"`), is done by computing expected value of the prediction over the distribution of random intercepts in parameter `logita`, or by using the two-point approximation described in Kangas et al. (XXXX).
`nrep`	The number of replicates in the Monte Carlo simulation when `varMethod="simul"`. Ignored if `varMethod="taylor"`.

Details

The variance-form-factor function is of form

$v(D,H,a,\lambda)=\pi \frac{\exp(a)}{1+\exp(a)} R(D,H,\lambda)^2 H$

where $a$ is the logit-transformed form factor and $R(D,H,\lambda)$ is the stem radius at stump height, which is approximated using

$R(D,H,\lambda)=w(H,\lambda)\frac{D}{2}+(1-w(H,\lambda))\frac{H}{H-B}\frac{D}{2}$

where the weight is taken from the right tail of the logit transformation

$w(H,\lambda)=2-2\frac{\exp\left( \frac{H-B}{\exp(\lambda)}\right)}{1+\exp\left( \frac{H-B}{\exp(\lambda)}\right)}$

Parameter uncertainty is reported because the same realized errors are used always when a model based on certain model fitting data is used; therefore those errors behave in practice like bias. Variance for total volume is computed as sum of all elements of the variance-covariance matrix of prediction errors of the mean.

Value

Function volvff returns a vector of tree volumes (in liters) that is of same length as vector dbh. In addition, attribute grad returns the Jacobian, which is used in nlme fitting for computing the derivatives of the model with respect to parameters, and in approximating the parameter uncertainty when varMethod="taylor".

Function predvff returns a list with following objects

`totvol`	Total volume of the trees of `dat`, in $m^3$ .
`totvolvar`	The estimated variance of `totvol`, taking into account parameter uncertainty.
`totvolci`	The estimated $100(1-p)$ % confidence interval of `totvol`, based on parameter uncertainty. If `varMethod="taylor"`, it is based on normality and the approximated variance `totvolvar`. If `varMethod="simul"`, the empirical quantiles of the `nrep` simulations are used.

And attributes

`trees`	A data frame of including tree-level volumes, their variance and 95% confidence intervals.
`varmu`	The variance-covariance matrix of prediction errors, taking into account theparameter uncertainty.

Author(s)

Lauri Mehtatalo <[email protected]>

References

Kangas A., Pitkanen T., Mehtatalo L., Heikkinen J. (xxxx) Mixed linear and non-linear tree volume models with regionally varying parameters

Examples

## Not run: 
library(lmfor)
data(treevol)

treevol$formfactor<-treevol$v/volvff(treevol$dbh,treevol$h,logita=100,lambda=log(0.2))
treevol$logitff<-log((treevol$formfactor)/(1-(treevol$formfactor)))

ptrees<-treevol[treevol$species=="pine",]
mod.init<-lm(logitff~I(1/h)+h+dbh+I(h*dbh)+I(1/(h*dbh))+
                     dataset+dataset:dbh+dataset:h,data=ptrees)

mod<-nlme(v~volvff(dbh,h,logita=logita,lambda=lambda),
         fixed=list(logita~I(1/h)+h+dbh+I(h*dbh)+I(1/(h*dbh))+dataset+
                           dataset:dbh+dataset:h+soil+temp_sum,
                           lambda~1),
         random=logita~1|stand/plot,
         start=c(coef(mod.init),rep(0,2),log(0.2)),
         data=ptrees,
         weights=varComb(varIdent(form=~1 |dataset),varPower()),
         method="ML",
#         control=list(msVerbose=TRUE),
#         verbose=TRUE
)

pred1<-predvff(ptrees,mod,varMethod="simul",biasCorr="integrate") 
pred1$totvol    
pred1$totvolvar
pred1$totvolci         
head(attributes(pred1)$trees)    

pred2<-predvff(ptrees,mod,varMethod="taylor",biasCorr="twopoint") 
pred2$totvol    
pred2$totvolvar
pred2$totvolci         
head(attributes(pred2)$trees)   

## End(Not run)
## Not run: 
library(lmfor)
data(treevol)

treevol$formfactor<-treevol$v/volvff(treevol$dbh,treevol$h,logita=100,lambda=log(0.2))
treevol$logitff<-log((treevol$formfactor)/(1-(treevol$formfactor)))

ptrees<-treevol[treevol$species=="pine",]
mod.init<-lm(logitff~I(1/h)+h+dbh+I(h*dbh)+I(1/(h*dbh))+
                     dataset+dataset:dbh+dataset:h,data=ptrees)

mod<-nlme(v~volvff(dbh,h,logita=logita,lambda=lambda),
         fixed=list(logita~I(1/h)+h+dbh+I(h*dbh)+I(1/(h*dbh))+dataset+
                           dataset:dbh+dataset:h+soil+temp_sum,
                           lambda~1),
         random=logita~1|stand/plot,
         start=c(coef(mod.init),rep(0,2),log(0.2)),
         data=ptrees,
         weights=varComb(varIdent(form=~1 |dataset),varPower()),
         method="ML",
#         control=list(msVerbose=TRUE),
#         verbose=TRUE
)

pred1<-predvff(ptrees,mod,varMethod="simul",biasCorr="integrate") 
pred1$totvol    
pred1$totvolvar
pred1$totvolci         
head(attributes(pred1)$trees)    

pred2<-predvff(ptrees,mod,varMethod="taylor",biasCorr="twopoint") 
pred2$totvol    
pred2$totvolvar
pred2$totvolci         
head(attributes(pred2)$trees)   

## End(Not run)

Package 'lmfor'

Help Index

Functions of Lauri Mehtatalo

Description

Details

Author(s)

References

Increment core data of Scots pine trees

Description

Usage

Format

Details

References

See Also

Examples

Individual tree characteristics and ALS data

Description

Usage

Format

Details

References

Examples

Breaking resistance (=bending strength) of birch wood samples

Description

Usage

Format

Details

References

Examples

Plot circles of a specified radius

Description

Usage

Arguments

Value

Author(s)

Examples

Evaluate the fit of a tree diameter distribution

Description

Usage

Arguments

Details

Value

Author(s)

Examples

Fit a Height-Diameter model to forest tree data using functions of package nlme.

Description

Usage

Arguments

Details

Value

Author(s)

References

See Also

Examples

CO2 exchange of transplanted Sphagnum fuscum moss in a chronosequence of mires.

Description

Usage

Format

Details

References

Examples

Available 2- and 3- parameter H-D model functions to be used by function fithd.

Description

Usage

Arguments

Details

Value

Author(s)

References

Examples

Estimate stand density using a Horvitz–Thompson-like estimator

Description

Usage

Arguments

Details

Value

Note

Author(s)

References

Examples

Fit a Height-Diameter model to forest tree data using functions of package `nlme`.

Available 2- and 3- parameter H-D model functions to be used by function `fithd`.