Package {morphErr}

morphErr: Measurement error models for morphometric data

Description

The morphErr package provides functions to analyse morphometric data using the method proposed by Stevenson, Smit, and Setyawan (2026).

Details

Two common ways of analysing morphometric data are to fit a linear regression model or find the reduced major axis (i.e., the principal component axis). However, these methods are known to perform poorly when observations are subject to non-negligible measurement error. The model implemented in this package provides the type of inference available from these methods, but explicitly accounts for measurement error. It can also handle multivariate data (rather than just fitting bivariate relationships), missing data, and accommodates dependence in measurement errors between dimensions.

Data requirements

The morphErr package might be right for you if your data were collected on photogrammetry surveys (e.g., using drones) and you have multiple photographs of some individuals.

The package might be particularly useful relative to alternative options if

• you take a lot of different measurements from each photograph (i.e., you observe a lot of "dimensions"),

• observations are subject to non-negligible measurement error,

• the measurement errors are correlated (e.g., a photograph with positive measurement error for one dimension tends to have positive measurement errors for other dimensions),

• you don't necessarily measure every dimension in every photograph,

• you'd like to estimate relationships amongst different subsets of the dimensions, or

• you'd like to predict different dimensions for different individuals using different subsets of the other dimensions.

Package overview

The key functions in morphErr are as follows:

• plotmorph() to plot morphometric data from a photogrammetry survey.

• fit.morph() to fit the model described by Stevenson, Smit, and Setyawan (2026).

• summary.lme.morph() with type = "pars" for parameter estimates and standard errors.

• summary.lme.morph() with type = "betas" for estimated coefficients of linear relationships to predict one dimension from any subset of the other dimensions.

• summary.lme.morph() with type = "betas-pca" for estimated coefficients of the reduced major axis (or pricipal component axis) summarising the relationship between two dimensions.

• summary.lme.morph() with type = "isometric-pca" or type = "isometric-pca-boot" to test for isometric growth between all pairs of dimensions.

• plot.lme.morph() to plot estimated relationships between dimensions.

• sim.measurements() to simulate morphometric data.

• sim.morph() to simulate multiple data sets and fit a model to each one.

• extract.sim.morph() to extract estimates from the models fitted using sim.morph().

Author(s)

Maintainer: Ben C. Stevenson ben.stevenson@st-andrews.ac.uk [copyright holder]

Authors:

• Ben C. Stevenson ben.stevenson@st-andrews.ac.uk [copyright holder]

• Elizabeth Smit esmi468@aucklanduni.ac.nz

References

Stevenson, B. C., Smit, E., and Setyawan, E. (2026) Measurement error models for morphometric data. Annals of Applied Statistics, 20: 945–962.

See Also

Useful links:

• https://github.com/elismit01/morphErr

• Report bugs at https://github.com/elismit01/morphErr/issues

Bootstrapping for a Morphometric Model

Description

Simulates bootstrap resamples and re-fits the original morphometric model to each one. The resulting object can be used as an argument for summary.lme.morph() and plot.lme.morph() to provide standard errors, confidence intervals, and p-values based on bootstrapping rather than asymptotic approximations.

Usage

boot.morph(object, B = 1000, control = NULL, progressbar = TRUE, n.cores = 1)

Arguments

Value

The original object, with additional bootstrapping information that can be used by summary.lme.morph() and plot.lme.morph()

Examples

## Fitting model to manta ray data.
fit <- fit.morph(manta)
## Carrying out a parametric bootstrap.
fit.boot <- boot.morph(fit, B = 100)
## Standard errors based on a bootstrap will be a little different
## to those based on asymptotic approximations.
summary(fit)
summary(fit.boot)

Calculate Beta Coefficients for Predictions

Description

Calculates coefficients for predicting one dimension from others, using either linear model or PCA interpretation.

Usage

calc.betas(
  fit,
  est = NULL,
  stders = TRUE,
  y.dim,
  x.dim,
  type = "lm",
  vcov = FALSE
)

Arguments

Value

If vcov = FALSE, a matrix with estimates and standard errors. If vcov = TRUE, a list with components:

Calculate Conditional Ratios Between Dimensions

Description

Calculates a ratio between dimensions, conditional on the denominator. Used internally by plotting functions.

Usage

calc.conditional.ratio(
  fit,
  stders = TRUE,
  y.dim,
  x.dim,
  newdata.x.dim,
  type = "lm"
)

Arguments

Details

Calculates ratios between dimensions conditional on the denominator dimension value. Can use either linear model or PCA interpretation for the calculations.

Value

A matrix with columns:

Construct Variance-Covariance Matrix

Description

Takes standard deviations and correlation parameters, and organises them into a block-diagonal variance-covariance matrix.

Usage

construct.varcov(sds, cors, n.blocks, m, block.only)

Arguments

Value

A variance-covariance matrix

Extract Results from a Morphometric Simulation Study

Description

Extracts results from each model fitted in a simulation study conducted using sim.morph().

Usage

extract.sim.morph(sim.res, FUN = summary, n.cores = 1)

Arguments

Value

An array containing the results of applying FUN to each model fit.

Examples

## Running a small simulation study.
sim.fits <- sim.morph(n.sims = 5,
                      n.animals = 10,
                      n.photos = 3,
                      mus = c(315, 150, 100),
                      sigmas = c(25, 15, 10),
                      rhos = c(0.85, 0.80, 0.75),
                      psis = c(10, 6, 4),
                      phis = c(0.5, 0.4, 0.3),
                      method = "REML",
                      progressbar = FALSE,
                      n.cores = 1)
## Extracting p-values from tests for isometry from each model
## fit.
iso.p <- function(x) summary(x, type = "isometric-pca")[, 3]
extract.sim.morph(sim.fits, FUN = iso.p)

Fit Morphometric Model

Description

Fits the model described by Stevenson, Smit, and Setyawan (2026) to morphometric data.

Usage

fit.morph(
  data,
  log.transform = FALSE,
  method = "REML",
  control = list(maxIter = 1e+05, msMaxIter = 1e+05)
)

Arguments

Details

This is a special case of a linear mixed-effects model, and is fitted via a call to nlme::lme(). Some arguments of fit.morph() are directly passed to nlme::lme().

Value

An object of classes lme.morph (specific to this package) and lme (inherited from nlme::lme()). The object is the list returned by nlme::lme() with a couple of additional components.

Rather than inspecting the object directly, the best way to extract understandable output from the object is using the S3 methods summary.lme.morph(), plot.lme.morph(), and predict.lme.morph() via the generic functions summary(), plot(), and predict(). Other S3 methods for the lme class are availble via nlme, such as nlme::coef.lme().

Log transformations

To fit a model to log-transformed measurement data, provide the untransformed measurements in the data argument and specify log.transform = TRUE. This will allow the S3 method plot.lme.morph() to plot the original data with back-transformed estimates, and the S3 method summary.lme.morph() to carry out the correct test for isometric growth.

The data argument

The manta object is an example of a correctly formatted data argument. It must be a data frame with the following columns:

animal.id

An individual identification number. Rows with the same animal.id correspond to measurements of the saime individual manta ray.

photo.id

A photo identification number. Rows with the same photo.id correspond to measurements taken from the same image.

photo.id

An integer indicating the dimension the measurement is for.

measurement

The corresponding measurement.

References

Stevenson, B. C., Smit, E., and Setyawan, E. (2026) Measurement error models for morphometric data. Annals of Applied Statistics, 20: 945–962.

Examples

## Fitting model to manta ray data.
fit <- fit.morph(manta)
## Using maximum likelihood instead of REML.
fit <- fit.morph(manta, method = "ML")

Manta Ray Morphometric Data

Description

Data collected on a drone photogrammetry survey of the reef manta ray Mobula alfredi in Raja Ampat, Indonesia.

Usage

manta

Format

manta

A data frame with four columns:

animal.id

An individual identification number. Rows with the same animal.id correspond to measurements of the same individual manta ray.

photo.id

A photo identification number. Rows with the same photo.id correspond to measurements taken from the same image.

dim

An integer indicating the dimension the measurement is for. Here, 1 is for disc width, 2 is for disc length, and 3 is for cranial width.

measurement

The observed measurement value in centimetres.

Organise Parameter Vector into Components

Description

Takes a vector of parameters and organises it into separate objects (e.g., a vector of mu parameters, variance-covariance matrices sigma and xi, etc).

Usage

organise.pars(pars, n.animals, n.photos, m, block.only = FALSE)

Arguments

Value

A list containing:

Plot Morphometric Data and Estimated Relationships

Description

An S3 method that plots morphometric data, estimated relationships between dimensions, or both, from a fitted model object returned by fit.morph().

Usage

## S3 method for class 'lme.morph'
plot(
  x,
  dims = c(1, 2),
  type = "data",
  line.type = "lm",
  confints = !add,
  add = FALSE,
  reverse.axes = FALSE,
  plot.data = TRUE,
  xlim = NULL,
  ylim = NULL,
  xlab = NULL,
  ylab = NULL,
  ...
)

Arguments

Value

No return value. Called to produce a plot.

Examples

## Fitting model to manta ray data.
fit <- fit.morph(manta)
## Plotting dimensions 1 and 2 with a fitted linear regression
## line.
plot(fit)
## Same again, but plotting the reduced major axis (or principal
## component axis) for dimensions 1 and 3.
plot(fit, dims = c(1, 3), line.type = "pca")
## A plot showing how the ratio between dimensions 2 and 3 changes
## with the size of dimension 3. Because the line is quite flat,
## it's plausible the relationship is isometric.
plot(fit, dims = c(3, 2), type = "ratio", line.type = "pca")
## On the other hand, the relationship between dimensions 1 and 2
## is clearly allometric.
plot(fit, dims = c(3, 1), type = "ratio", line.type = "pca")

Plot Morphometric Data

Description

Creates scatter plots of morphometric measurements, with options for different dimension combinations and plotting of ratios.

Usage

plotmorph(
  data,
  dims = c(1, 2),
  plot.data = TRUE,
  ratios = FALSE,
  xlim = NULL,
  ylim = NULL,
  xlab = NULL,
  ylab = NULL
)

Arguments

Value

No return value. Called to produce a plot.

The data argument

The manta object is an example of a correctly formatted data argument. It must be a data frame with the following columns:

animal.id

An individual identification number. Rows with the same animal.id correspond to measurements of the saime individual manta ray.

photo.id

A photo identification number. Rows with the same photo.id correspond to measurements taken from the same image.

photo.id

An integer indicating the dimension the measurement is for.

measurement

The corresponding measurement.

Examples

## Plotting dimensions 1 and 2.
plotmorph(manta)
## Plotting dimensions 1 and 3.
plotmorph(manta, dims = c(1, 3))
## Plotting the ratio of dimension 2 divided by dimension 1 on the
## y-axis.
plotmorph(manta, dims = c(1, 2), ratios = TRUE)

Predict Measurements from True Values

Description

Calculates model predictions for the true size of a dimension based on true values for any subset of the remaining dimensions. Predictions for multiple individuals are available from a single call to the function.

Usage

## S3 method for class 'lme.morph'
predict(object, y.dim, newdata = NULL, type = c("lm", "pca"), ...)

Arguments

Details

This function computes estimated coefficients for the appropriate linear function (as per summary.lme.morph() with type = "betas-lm" or type = "betas-pca"), and then evaluates the function for the provided newdata. For type = "pca", only a single column can be provided in newdata.

Value

A matrix with estimated true sizes and standard errors. Note that these are standard errors, not prediction errors.

Prediction functions in morphErr

There are three key differences between predict.lme.morph() and predictblup().

• predict.lme.morph() only generates estimates for one dimension based on true values for other dimensions, whereas predictblup() can generate estimates from true values, observed values (i.e., subject to measurement error), or a combination of both.

• predict.lme.morph() can generate estimates for multiple individuals, whereas predictblup() only provides estimates for a single individual.

• predict.lme.morph() provides standard errors, but predictblup() does not.

See Also

predictblup()

Examples

## Fitting model to manta ray data.
fit <- fit.morph(manta)
## Predicting dimension 2 for a single individual with a true value
## of 90 for dimension 3.
predict(fit, y.dim = 2, newdata = data.frame(dim3 = 90))
## Predicting dimension 1 from dimensions 2 and 3 for two
## individuals, one with true values of 130 and 60 for dimensions 2
## and 3, respectively, and one with true values of 140 and 70.
predict(fit, y.dim = 1, newdata = data.frame(dim2 = c(130, 140),
                                             dim3 = c(60, 70)))

Predict from Observed Measurements

Description

Makes predictions for a single individual using the best linear unbiased predictor (BLUP).

Usage

predictblup(object, true = NULL, obs = NULL)

Arguments

Details

This function uses a BLUP to compute estimated true values for a single individual. A BLUP is computed by finding the mode of the multivariate probability density function of the true values, conditional on any provided values for true dimension sizes, dimension measurements observed with error, or a combination of both.

Value

A numeric vector of predictions for all dimensions.

Prediction functions in morphErr

There are three key differences between predict.lme.morph() and predictblup().

• predict.lme.morph() only generates estimates for one dimension based on true values for other dimensions, whereas predictblup() can generate estimates from true values, observed values (i.e., subject to measurement error), or a combination of both.

• predict.lme.morph() can generate estimates for multiple individuals, whereas predictblup() only provides estimates for a single individual.

• predict.lme.morph() provides standard errors, but predictblup() does not.

See Also

predict.lme.morph()

Examples

## Fitting model to manta ray data.
fit <- fit.morph(manta)
## Estimates for the true dimension sizes of an individual with a
## known true value for dimension 2 of 130, and observed 
## measurements subject to error from two photographs. The first 
## photograph has # observed measurements of 300 and 135 for 
## dimensions 1 and 2, respectively, while the second photograph 
## has observed # measurements of 290 and 140, respectively.
## Dimension 3 is not observed in any photograph.
predictblup(fit, true = c(NA, 130, NA), obs = rbind(c(300, 135, NA),
                                                    c(290, 140, NA)))

Simulate Morphometric Measurements

Description

Simulates morphometric data from a photogrammetry survey, with observations subject to measurement error, under the model described by Stevenson, Smit, and Setyawan (2026).

Usage

sim.measurements(
  n.animals = NULL,
  n.photos = NULL,
  data = NULL,
  mus,
  sigmas,
  rhos,
  psis,
  phis,
  log.transform = FALSE
)

Arguments

Details

For arguments rhos and phis, the elements must be ordered so that all m - 1 correlations involving dimension 1 appear first in ascending numerical order, followed by all remaining m - 2 correlations involving dimension 2, and so on, where ⁠eqn{m}⁠ is the number of dimensions.

For example, if m = 4, then the first three elements are the correlations between dimension 1 and dimensions 2, 3, and 4, respectively. The following two elements are correlations between dimension 2 and dimensions 3 and 4, respectively. The final element is the correlation between dimension 3 and 4.

Value

A data frame with four columns:

animal.id

An individual identification number. Rows with the same animal.id correspond to measurements of the same individual.

photo.id

A photo identification number. Rows with the same photo.id correspond to measurements taken from the same image.

dim

An integer indicating the dimension the measurement is for.

measurement

The observed measurement value.

References

Stevenson, B. C., Smit, E., and Setyawan, E. (2026) Measurement error models for morphometric data. Annals of Applied Statistics, 20: 945–962.

See Also

sim.morph() to conduct a simulation study by simulating multiple data sets and fitting a model to each one.

Examples

## Simulating data for ten animals, with two photos each, measuring
## three dimensions.
sim.data <- sim.measurements(n.animals = 10, n.photos = 2,
                             mus = c(315, 150, 100),
                             sigmas = c(25, 15, 10),
                             rhos = c(0.85, 0.80, 0.75),
                             psis = c(10, 6, 4),
                             phis = c(0.5, 0.4, 0.3))
head(sim.data)

## Simulating data for two animals, with different numbers of
## photos for each, and different measurements available from
## different photos.
data <- data.frame(animal.id = c(1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2),
                   photo.id = c(1, 1, 2, 2, 2, 1, 2, 3, 3, 4, 4),
                   dim = c(1, 3, 1, 2, 3, 1, 1, 1, 3, 2, 3))
sim.data <- sim.measurements(data = data,
                             mus = c(315, 150, 100),
                             sigmas = c(25, 15, 10),
                             rhos = c(0.85, 0.80, 0.75),
                             psis = c(10, 6, 4),
                             phis = c(0.5, 0.4, 0.3))
sim.data

Conduct a Simulation Study for Morphometric Models

Description

Simulates multiple datasets and fits a model to each one.

Usage

sim.morph(
  n.sims,
  n.animals = NULL,
  n.photos = NULL,
  data = NULL,
  mus,
  sigmas,
  rhos,
  psis,
  phis,
  log.transform = FALSE,
  method = "REML",
  control = list(maxIter = 1e+05, msMaxIter = 1e+05),
  progressbar = TRUE,
  n.cores = 1
)

Arguments

Value

An object of class lme.morph.sim. The best way to extract results from this object is to use extract.sim.morph(). See the example below.

See Also

extract.sim.morph() to extract results from the object returned by this function.

Examples

## Running a small simulation study.
sim.fits <- sim.morph(n.sims = 5,
                      n.animals = 10,
                      n.photos = 3,
                      mus = c(315, 150, 100),
                      sigmas = c(25, 15, 10),
                      rhos = c(0.85, 0.80, 0.75),
                      psis = c(10, 6, 4),
                      phis = c(0.5, 0.4, 0.3),
                      method = "REML",
                      progressbar = FALSE,
                      n.cores = 1)
## Extracting p-values from tests for isometry from each model
## fit.
iso.p <- function(x) summary(x, type = "isometric-pca")[, 3]
extract.sim.morph(sim.fits, FUN = iso.p)

Summarise Morphometric Model Fits

Description

An S3 method that creates summaries of fitted morphometric models. This function can provide different types of summaries including parameter estimates, beta coefficients, and tests for isometry.

Usage

## S3 method for class 'lme.morph'
summary(object, ..., type = "pars", y.dim, x.dim, level = 0.95)

Arguments

Details

The following summaries are available with this method, controlled by the type argument:

type = "pars"

Estimates and standard errors for the model parameters, including the vectors \mu for the means of the true dimension sizes, \sigma for the standard deviations for true dimension sizes, \rho for the correlations between true dimension sizes, \psi for the standard deviations in measurement errors for the dimensions, and \phi for the correlations between measurement errors for all pairs of dimensions.

type = "betas" or type = "betas-lm"

Estimates and standard errors for coefficients of a linear combination that provides the expected value of one dimension (specified by y.dim) using any subset of the remaining dimensions (specified in the vector x.dim).

type = "betas-pca"

An estimated intercept and slope for the reduced major axis (or principal component axis) summarising the relationship between the true values of dimensions specified by y.dim and x.dim.

type = "isometric-pca"

Results for tests of the null hypotheses of isometric relationships between all pairs of dimensions. The p-values are obtained using normal approximations for the test statistics.

For models with log.transform = FALSE, the test statistic to test for isometric growth between dimensions p and q is \mu_p\sigma_q - \mu_q\sigma_p, which is equal to 0 under the null hypothesis.

For models with log.transform = TRUE, the test statistic is the slope for the reduced major axis (or principal component axis) summarising the relationship between the two dimensions, which is equal to 1 under the null hypothesis.

Value

A matrix or data frame containing the requested summary.

Examples

## Fitting model to manta ray data.
fit <- fit.morph(manta)
## Parameter estimates and standard errors.
summary(fit)
## Estimated coefficients for the linear combination that provides
## the expected value of the first dimension from only the second.
summary(fit, type = "betas-lm", y.dim = 1, x.dim = 2)
## Estimated coefficients for the linear combination that provides
## the expected value of the second dimension from both the first
## and third.
summary(fit, type = "betas-lm", y.dim = 2, x.dim = c(1, 3))
## Estimated intercept and slope for the reduced major axis (or
## pricipal component axis) summarising the relationship between
## the first and second dimensions.
summary(fit, type = "betas-pca", y.dim = 1, x.dim = 2)
## Tests for isometry between all pairs of dimensions.
summary(fit, type = "isometric-pca")

Extract Variance-Covariance Matrix from Morphometric Model

Description

An S3 method for extracting both parameter estimates and their variance-covariance matrix from a fitted morphometric model.

Usage

## S3 method for class 'lme.morph'
vcov(object, ...)

Arguments

Value

A list containing:

Examples

## Fitting model to manta ray data.
fit <- fit.morph(manta)
## Extracting variance-covariance matrix.
vcov(fit)