Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.
In a recent post I have shown that we can build linear combinations of model parameters (see here ). For example, if we have two parameter estimates, say
where
In biology, nonlinear transformations are much more frequent than linear transformations. Nonlinear transformations are, e.g.,
Let’s start with the transformation of a single variable
Example 1: The half-life of a herbicide
A herbicide has proven to follow a first order degradation kinetic in soil, with constant degradation rate
Every pesticide chemist knows that the half-life is derived by the degradation rate, according to the following equation:
Therefore, if we plug-in x = X, the half-life of our herbicide is:
Y <- log(0.5)/-0.035 Y
[1] 19.80421
But … what is the standard error of this half-life? There is some uncertainty around the estimate of
The basic idea beyond the delta method
The basic idea behind the delta method is that most of the simple nonlinear functions, which we use in biology, can be locally approximated by the tangent line through a point of interest. For example, our nonlinear half-life function is
What is the equation of the tangent line? In general, if the nonlinear function is
We need the derivative!
In order to write the equation of the red line in the Figure above, we need to consider that D(). For our case, it is:
D(expression(log(0.5)/X), "X")
-(log(0.5)/X^2)
Therefore, we can use this R function to calculate the slope
X <- -0.035 m <- eval( D(expression(log(0.5)/X), "X") ) m
[1] 565.8344
We already know that
that is:
Replacing a curve with a line
Now, we have two functions:
- the original nonlinear half-life function
$ - a new linear function (
), that is a very close approximation to the previous one, at least near to the point , which we are interested in.
Therefore, we can approximate the former with the latter! If we use the linear function, we see that the half-life is:
39.60841 + 565.8344 * -0.035
[1] 19.80421
which is what we expected. The advantage is that we can now use the low of propagation of errors to estimate the standard error (see the first and second equation in this post):
Here we go:
sqrt( m^2 * (0.00195 ^ 2) )
[1] 1.103377
In general…
If we have a nonlinear transformation
Example 2: a back-transformed count
A paper reports that the mean number of microorganisms in a substrate, on a logarithmic scale, was
The first derivative of our nonlinear function is:
D(expression(exp(X)), "X")
exp(X)
and thus the slope of the tangent line is:
X <- 5 m <- eval( D(expression(exp(X)), "X") ) m
[1] 148.4132
According to the function above, the standard error for the back-transformed mean is:
sigma <- 0.84 sqrt( m^2 * sigma^2 )
[1] 124.6671
Example 3: Selenium concentration in olive drupes
When the combination is nonlinear and involves two parameters, such as:
the delta method requires a slightly more complex formula, which we will see by using another example.
The concentration of selenium in olive drupes was found to be
The amount of selenium is easily calculated as:
X <- 3.1; Z = 3.4 Y <- X * Z Y
[1] 10.54
Delta standard errors can be obtained by considering the partial derivatives for each of the two variables and combining them as follows:
dX <- eval( D(expression(X * Z), "X") ) dZ <- eval( D(expression(X * Z), "Z") ) sigmaX <- 0.8; sigmaZ <- 0.31; sigmaXZ <- 0.55 sqrt( (dX^2) * sigmaX^2 + (dZ^2) * sigmaZ^2 + 2 * dX * dZ * sigmaXZ )
[1] 4.462726
In general terms, when
For those of you who would like to get involved with matrix notation: we can reach the same result via matrix multiplication (see below). This might be easier when we have more than two variables to combine.
der <- matrix(c(dX, dZ), 1, 2) sigma <- matrix(c(sigmaX^2, sigmaXZ, sigmaXZ, sigmaZ^2), 2, 2, byrow = T) sqrt( der %*% sigma %*% t(der) )
[,1] [1,] 4.462726
The delta method with R
In R there is a shortcut function to calculate delta standard errors, that is available in the ‘car’ package. In order to use it, we need to have:
- a named vector for the variables that we have to combine
- an expression for the transformation
- a variance-covariance matrix
For the first example, we have:
obj <- c("X" = -0.035)
sigma <- matrix(c(0.00195^2), 1, 1)
library(car)
Loading required package: carData
deltaMethod(object = obj, g="log(0.5)/X", vcov = sigma)
Estimate SE 2.5 % 97.5 % log(0.5)/X 19.8042 1.1034 17.6416 21.967
For the second example:
obj <- c("X" = 5)
sigma <- matrix(c(0.84^2), 1, 1)
deltaMethod(object = obj, g="exp(X)", vcov = sigma)
Estimate SE 2.5 % 97.5 % exp(X) 148.41 124.67 -95.93 392.76
For the third example:
obj <- c("X" = 3.1, "Z" = 3.4)
sigma <- matrix(c(0.8^2, 0.55, 0.55, 0.31^2), 2, 2, byrow = T)
deltaMethod(object = obj, g="X * Z", vcov = sigma)
Estimate SE 2.5 % 97.5 % X * Z 10.5400 4.4627 1.7932 19.287
The function ‘deltaMethod()’ is very handy to be used in connection with model objects, as we do not need to provide anything, but the transformation function. But this is something that requires another post!
However, two final notes relating to the delta method need to be pointed out here:
- the delta standard error is always approximate;
- if the original variables are gaussian, the transformed variable, usually, is not gaussian.
Thanks for reading—and don’t forget to check out my new book below!
Andrea Onofri
Department of Agricultural, Food and Environmental Sciences
University of Perugia (Italy)
Send comments to: andrea.onofri@unipg.it
Reference
- Beyene, J., Moineddin, R., 2005. Methods for confidence interval estimation of a ratio parameter with application to location quotients. BMC Med Res Methodol 5, 32. https://doi.org/10.1186/1471-2288-5-32
Post written on 25/5/2019 and edited on 22/11/2024 and 12/12/2025
R-bloggers.com offers daily e-mail updates about R news and tutorials about learning R and many other topics. Click here if you're looking to post or find an R/data-science job.
Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.
