## Friday, August 17, 2012

### The DELTA method

# The delta method is an extremely useful tool for estimating the standard errors of non-linear models.

# This post will use as a reference David Patterson of the University of Montana's post on the delta method http://www.math.umt.edu/patterson/549/Delta.pdf

# The delta method states that the standard error of a distribution f(x) is approximated by f'(Mx) var(x) f(Mx) where f(x) is the gradiant of F(x) with respect to x (first derivative), and Mx=mean(x).

# Let's first generate a base variable x~N with mean 5 and standard deviation 3.
x=rnorm(10000)*3+5

# Example 1: F(x) = x^2 -> f(x) = 2x ->
# thus var(F(x))~2*mean(x)*var(x)*2*mean(x) ->
# var(F(x)) ~
2*5*9*2*5

# Let's see how well our approximation compares with actual data
x2 = x^2
var(x2)

# We can see that our estimate (900) is pretty close to the sample variance (1050ish)

# Using our sampling distribution
2*mean(x)*var(x)*2*mean(x)

# Example 2: F(x) = 1/x -> f(x) = (-1)*x^(-2)
# var(F(x)) ~ (-1)*mean(x)^(-2) * var(x) * (-1)*mean(x)^(-2)
# = mean(x)^(-4) * var(x)
mean(x)^(-4) * var(x)

# Let's see how well our approximation compares with actual data
xi = 1/x

var(xi)
# In this case, not very well.  This is probably because the distribution of x approaches the nondiferentiable point 0 for some xs making a first order approximation fail badly .

# Let us try with another variable
z = rnorm(10000)+20

mean(z)^(-4) * var(z)

zi = 1/z

var(zi)
# We can see now the delta method is working well.

# Example 3: F(z) = ln(z) -> f(z) = 1/z
# var(F(z)) ~ (1/mean(z))*var(z)*(1/mean(z))
# = (1/mean(z))^2*var(z)
(1/mean(z))^2*var(z)

lnz = log(z)

var(lnz)
# Once again, looking good.

# Now let's imagine an estimator.  First off we need to recognize that estimators are random variables.

# Bhat = (X'X)^(-1)X'Y = (X'X)^(-1)X'(XB + e) = (X'X)^(-1)X'XB + (X'X)^(-1)X'e
#      = B + (X'X)^(-1)X'e
# E(Bhat|x) = E(B|x) + E[(X'X)^-1 X'e|x] -> assume (e|x)=0 ->
# E(Bhat|x) = B + (X'X)^-1 X'E[e|x] = B
# var(Bhat|x) = var(B + (X'X)^(-1)X'e|x) = var((X'X)^(-1)X'e|x)
# = (X'X)^(-1)X'var(e|x)((X'X)^(-1)X')' = (X'X)^(-1)X'var(e|x)X(X'X)^(-1)
# = (X'X)^(-1)X' *I X(X'X)^(-1)
# = sigma2*(X'X)^(-1)X'X(X'X)^(-1) = sigma2*(X'X)^(-1)

e = rnorm(1000)
y = 5*x + e*100

# Let's estimate OLS without residuals
ols=lm(y~x-1)

summary(ols)

uhat = residuals(ols)

# var(Bhat|x) ~ sigma2(X'X)^(-1) estimate
var(uhat)*sum(x*x)^(-1)

# standard error
(var(uhat)*sum(x*x)^(-1))^.5

#### 1 comment:

1. Good example of being careful when using the delta method. A note to readers that are wondering if you always have to do the delta method manually: the car package has a nice function 'deltaMethod'.