# Custom power analysis in R

Distance sampling is a common technique used to estimate the abundance of wildlife species (such as deer or turkey). This technique employees line-transect based sampling where perpendicular distance to individual animals is measured. Those distances are then used to model a detection function. That detection function is then used to correct counts of animals for incomplete detection (Buckland et al. 2001). Wildlife surveys that employee distance sampling used Program Distance (Thomas et al. 2010) to estimate abundance. Program Distance results in an abundance estimate $\hat N$ and it's coefficient of variation $cv(\hat N)$ where $cv(\hat N)$ is $se(\hat N)/\hat N$ and $var(\hat N)$ is $se(\hat N)^2$. To compare abundances from two different areas, one can use a simple z-test: $$Z = \frac{(\hat N_1-\hat N_2)}{\sqrt{(var(\hat N_1-\hat N_2))}}$$ where $$var(\hat N_1-\hat N_2) = var(\hat N_1)+var(\hat N_2)$$ provided $\hat N_1$ and $\hat N_2$ were independently estimated (Buckland et al. 2001).

After conducting a pilot survey, it is easy to estimate the amount of effort (i.e., length of transects $\hat L$) needed to obtain a target coefficient of variation desired from future survey efforts. $$\hat L = \frac{(L_o*cv(\hat N)^2)}{(cv_t(\hat N)^2)}$$ where $L_o$ is the total line length surveyed during a pilot survey, $cv(\hat N)$ is the coefficient of variation of the abundance estimate from the pilot survey, and $cv_t(\hat N)$ is the target coefficient of variation desired from future survey efforts (Buckland et al. 2001).

To determine a target coefficient of variation, I have developed a R function to estimate the power of the above test to detect a given difference between $\hat N_1$ and $\hat N_2$.

# A function to estimate power to detect a change between two abundance estimates
# using a 2-tailed z-test. Where N1 and N2 are the two abundances, change is the
# desired percent change or difference, and CV is the coefficient of variation.
ZPOWER.2tail = function(change,N1,CV,alpha){
N2=N1+(N1*change)
SD1=N1*CV
SD2=N2*CV
VAR=(SD1^2)+(SD2^2)
qL=qnorm(alpha/2)
qU=qnorm(1-(alpha/2))
Z=(N1-N2)/sqrt(VAR)
power=(pnorm(qL-Z))+(1-pnorm(qU-Z))
cbind(N1,N2,change,CV,SD1,SD2,VAR,qL,qU,Z,alpha,power)}

# Example use of function.
ZPOWER.2tail(change=-0.2,N1=300,CV=0.15,alpha=0.05)


A one-tailed verison would be:

# A function to estimate power to detect a change between two abundance estimates
# using a 1-tailed z-test. Where N1 and N2 are the two abundances, change is the
# desired percent change or difference, and CV is the coefficient of variation.
ZPOWER.1tail = function(change,N1,CV,alpha){
N2=N1+(N1*change)
SD1=N1*CV
SD2=N2*CV
VAR=(SD1^2)+(SD2^2)
if(change>0) q=qnorm(alpha)
if(change<0) q=qnorm(1-alpha)
z=(N1-N2)/sqrt(VAR)
if(change>0) power=pnorm(q-z)
if(change<0) power=1-pnorm(q-z)
cbind(N1,N2,change,CV,SD1,SD2,VAR,q,z,alpha,power)}


Of course these functions could easily be combined in to one function after all the bugs are worked out.

After writing this code, I found the following code posted at CREEM.

#--------------------------------------------
# Power to detect a change between two time
# points using a z-test.
#--------------------------------------------
#alpha level (type I error)
alpha<-0.05
#vector of coefficients of var
cvs <- seq(0.05,0.5,by=0.05)
#vector of true differences
fs <- seq(0,1,by=0.1)
lcvs<-length(cvs)
lfs<-length(fs)
res<-matrix(0,lcvs,lfs)
q1<-qnorm(1-(alpha/2))
q2<-qnorm(alpha/2)
s2<-sqrt(2)
#work out power at each level of CV and true diff
for(i1 in 1:lcvs) {
cv<-cvs[i1]
for (i2 in 1:lfs) {
f<-fs[i2]
z<-f/(s2*cv)
power = (1-pnorm(q1-z)) + pnorm(q2 - z)
res[i1,i2]<-power
}
}
#output as a perspective plot
persp(cvs,fs,res)
#--------------------------------------------


I know it was a long time getting here, but here is my question.

The major difference between my code and the CREEM code is: $$Z=(\hat N_1-\hat N_2)/\sqrt{var(\hat N_1-\hat N_2}$$ vs. $$Z=d/(\sqrt{2}*cv(\hat N))$$ Where $d$ is the "true difference."

Can someone help me understand the use of $\sqrt{2}*cv(\hat N)$ in the denominator instead of $\sqrt{var(\hat N_1-\hat N_2)}$ ?

Thanks.

Buckland, S. T., D. R. Anderson, K. P. Burnham, J. L. Laake, D. L. Borchers, and L. Thomas. 2001. Introduction to distance sampling: estimating abundance of biological populations. Oxford University Press, New York, New York, USA.

Thomas, L., S. T. Buckland, E. A. Rexstad, J. L. Laake, S. Strindberg, S. L. Hedley, J. R. B. Bishop, T. A. Marques, and K. P. Burnham. 2010. Distance software: design and analysis of distance sampling surveys for estimating population size. Journal of Applied Ecology 47:5–14.

-

If one assumes $d = \frac{(\hat N_2-\hat N_1)}{\hat N_2}$, which is the percent differnce. Then:
$$Z = \frac{d}{\sqrt{2}*cv(\hat N)}$$ $$Z = \frac{\frac{(\hat N_2-\hat N_1)}{\hat N_2}}{\sqrt{2}*cv(\hat N)}$$ $$Z = \frac{(\hat N_2-\hat N_1)}{\sqrt{2}*cv(\hat N)*\hat N_2}$$ $$Z = \frac{(\hat N_2-\hat N_1)}{\sqrt{2}*se(\hat N)}$$ $$Z = \frac{(\hat N_2-\hat N_1)}{\sqrt{se(\hat N)^2+se(\hat N)^2}}$$ Thus, the CREEM function assumes that $se(\hat N_1)$ = $se(\hat N_2)$ but my fuction assumes $cv(\hat N_1)$ = $cv(\hat N_2)$.
Therefore, the primary difference is in the assumption about the relationship between variance and abundance. For example, one could assume variance to be proportional to $\hat N$, $\hat N^2$, or $\hat N^3$ (Gerrodette 1987) or constant as in the case of the CREEM function. Gerrodette (1987) suggested that constant $cv(\hat N)$ was an appropriate assumption for distance sampling-based estimates of abundance. Estimates from mark-recapture based surveys might be better suited to the assumption that $cv(\hat N)$ is proportional to $\sqrt {\hat N}$.