sma: (Standardized) major axis estimation and testing for one or several samples

Description

The sma and ma functions fit SMA and MA respectively, and construct confidence intervals or test hypotheses about slope or elevation in one or several samples, depending on how the arguments are specified. Options exist to force lines through the origin, (approximately) correct for measurement error if the measurement error variance is known, and use robust estimation to ensure that inferences are valid in the presence of outliers (currently implemented for single samples only).

Usage

sma(formula, data, subset, na.action, log="",
	 method=c("SMA","MA","OLS"), type=c("elevation","shift"), 
   alpha=0.05, slope.test=NA, elev.test=NA, 
   multcomp=FALSE, multcompmethod=c("default","adjusted"),
	 robust=FALSE,V=matrix(0,2,2),n_min=3,quiet=FALSE,
	 ...)

Arguments

formula

Formula of the form y ~ x etc. This determines whether a single (S)MA is fitted, or whether several lines are fitted and compared. See Details.

data

A dataframe with the x and y variables.

subset

A subset of the dataframe for fitting; optional.

na.action

What to do with missing values (na.omit, na.fail, etc).

log

One of 'x', 'y', or 'xy' to log10-transform variables.

method

If SMA, standardized major axis, if MA, major axis, and if OLS, ordinary least squares.

type

If several lines with common slope are to be compared, do you want to test for a change in 'elevation' or for a 'shift' along a common (S)MA. See Details.

alpha

The error rate for confidence intervals. Typically 0.05.

slope.test

The hypothesised value to be used, if testing for evidence that (S)MA slope(s) are significantly different from a hypothesised value.

elev.test

The hypothesised value to be used, if testing for evidence that (S)MA elevation(s) are significantly different from a hypothesised value.

robust

If TRUE, uses a new method of robust line fitting. (Currently available for single-sample testing only.)

The estimated variance matrix of measurement error. Average measurement error for Y is in the first row and column, and average measurement error for X is in the second row and column. The default is that there is no measurement error.

n_min

The minimum sample size for a group.

quiet

If TRUE, suppresses all messages.

multcomp

Logical. If TRUE, performs pair-wise comparisons between levels of the grouping variable.

multcompmethod

Whether to adjust the P values for multiple comparisons ('adjusted') or not ('default').

…

Further arguments passed to internal functions (none at the moment?)

Value

An object of class sma or ma, which contains the following output arguments:

coef: The coefficients of the fitted (standardised) major axes. If several samples are being compared, this will return parameters from the alternative model, e.g. assuming separate slopes if testing for common slope, or assuming common slope but separate elevations if testing for common elevation.
nullcoef: The coefficients under the null hypothesis
alpha: As above.
method: The method used in fitting lines: 'MA' or 'SMA'
intercept: Whether or not (S)MA lines were forced through the origin (True or False).
call: The call to the ma or sma function.
data: As above.
log: As above.
variables: A list of the variables used in fitting (S)MA lines.
origvariables: A list of the variables prior to transformation (if any) for use in fitting (S)MA lines.
groups: Levels of the grouping variable, if present in the fit.
gt: Type of grouptest ("slopecom","elevcom",or "shiftcom"), if it was carried out, or "none" if none.
gtr: The result of that grouptest.
slopetest: Output from the hypothesis test of slope(s), if any. Returned as a list of objects, including the P-value (p), the test statistic (r or LR), the (common) slope (b) and its confidence interval (ci).
elevtest: Output from the hypothesis test of elevation(s), if any. Returned as a list of objects, including the P-value (p), the test statistic (t or stat), the (common) elevation (a) and its confidence interval (ci).
slopetestdone: Whether a slopetest was actually carried out.
elevtestdone: Whether an elevation test was actually carried out.
n: Sample size(s).
r2: Squared correlation coefficient.
pval: P-value of the test of correlation coefficient against zero.
from,to: X values corresponding the maximum and minimum fitted values in each group. Used by plot.sma to determine endpoints for fitted lines).
groupsummary: Neatly organized dataframe with coefficients by group.

Details

This is the key function in the smatr-package; all the key estimation and testing functionality in the package can all be accessed using this function, via different usages of the formula and other arguments, as described below.

One-sample testing The below options allow estimation of a (S)MA, confidence intervals for parameters, and hypothesis testing of parameters, from a single sample of two variables y and x. Use the sma function to fit a standardised major axis (SMA), or use ma in combination with the below options in order to fit major axis (MA) instead.

sma(y~x): Fits a SMA and constructs confidence intervals for the true slope and elevation. Replaces the line.cis function from previous versions of smatr.
ma(y~x): Fits a MA and constructs confidence intervals for the true slope and elevation. All the below functions also work for MA, if the ma function is called instead of the sma function.
sma(y~x, slope.test=B): Tests if the slope of a SMA equals B.
sma(y~x, elev.test=A): Tests if the elevation of a SMA equals A.
sma(y~x, robust=T): Fits a SMA using Huber's M estimation and constructs confidence intervals for the true slope and elevation. This offers robustness to outliers in estimation and inference, and can be used in combination with the slope.test and elev.test arguments.
sma(y~x-1): Fits a SMA where the line is forced through the origin, and constructs confidence intervals for the true slope. This type of formula can be used in combination with the slope.test argument.

For several samples: The below options allow estimation of several (S)MA lines, confidence intervals for parameters, and hypothesis testing of parameters, from two variables y and x for observations that have been classified into several different samples using the factor groups. Use the sma function to fit a standardised major axis (SMA), or use the ma in combination with the below options in order to fir major axis (MA) instead.

sma(y~x*groups): Test if several SMA lines share a common slope, and construct a confidence interval for the true common slope.
sma(y~x+groups, type="elevation"): Test if several common slope SMA lines also share a common elevation, and construct a confidence interval for the true common elevation.
sma(y~x+groups, type="shift"): Test if several groups of observations have no shift in location along common slope SMA lines.
sma(y~x*groups, slope.test=B): Test if several SMA lines share a common slope whose true value is exactly equal to B.
sma(y~x+groups, elev.test=A): Test if several common-slope SMA lines share a common elevation whose true value is exactly equal to A.
sma(y~x*groups-1): Test if several SMA lines forced through the origin share a common slope. This can also be used in combination with the slope.test argument or when testing for no shift along common (S)MA lines.

In all cases, estimates and confidence intervals for key parameters are returned, and if a hypothesis test is done, results will be returned and stored in the slope.test or elev.test output arguments.

The plot function can be applied to objects produced using the sma and ma functions, which is highly recommended to visualise results and check assumptions.

Multiple comparisons If multcomp=TRUE, pair-wise comparisons are made between levels of the grouping variable for differences in slope, elevation or shift, depending on the formula used. The P values can be adjusted for multiple comparisons (using the Sidak correction). See also multcompmatrix for visualization of the results. Warning: When using the multiple comparisons (multcomp=TRUE), you must specify a data statement. If your variables are not in a dataframe, simply combine them in a dataframe before calling sma.

References

Warton, D.I., R.A. Duursma, D.S. Falster and S. Taskinen (2012). smatr 3 - an R package for estimation and inference about allometric lines. Methods in Ecology and Evolution. 3, 257-259.

Warton D. I. and Weber N. C. (2002) Common slope tests for bivariate structural relationships. Biometrical Journal 44, 161--174.

Warton D. I., Wright I. J., Falster D. S. and Westoby M. (2006) A review of bivariate line-fitting methods for allometry. Biological Reviews 81, 259--291.

Taskinen S. and Warton D. I. (in press) Robust estimation and inference for bivariate line-fitting in allometry. Biometrical Journal.

Examples

Run this code

# NOT RUN {
# Load leaf lifetime dataset:
data(leaflife)

### One sample analyses ###
# Extract only low-nutrient, low-rainfall data:
leaf.low <- subset(leaflife, soilp == 'low' & rain == 'low')

# Fit a MA for log(leaf longevity) vs log(leaf mass per area):
ma(longev ~ lma, log='xy', data=leaflife)

# Test if the MA slope is not significantly different from 1:
ma.test <- ma(longev ~ lma, log='xy', slope.test=1, data=leaflife)
summary(ma.test)

# Construct a residual plot to check assumptions:
plot(ma.test,type="residual")

### Several sample analyses ###

# Now consider low-nutrient sites (high and low rainfall):
leaf.low.soilp <- subset(leaflife, soilp == 'low')

# Fit SMA's separately at each of high and low rainfall sites,
# and test for common slope:
com.test <- sma(longev~lma*rain, log="xy", data=leaf.low.soilp)
com.test

# Plot longevity vs LMA separately for each group:
plot(com.test)

# Fit SMA's separately at each of high and low rainfall sites,
# and test if there is a common slope equal to 1:
sma(longev~lma*rain, log="xy", slope.test=1, data=leaf.low.soilp)

# Fit SMA's with common slope across each of high and low rainfall sites, 
# and test for common elevation:
sma(longev~lma+rain, log="xy", data=leaf.low.soilp)

# Fit SMA's with common slope across each of high and low rainfall sites, 
# and test for no shift along common SMA:
sma(longev~lma+rain, log="xy", type="shift", data=leaf.low.soilp)

# }

Run the code above in your browser using DataLab