Variance of a parameter average on the real scale

posts related to the RMark library, which may not be of general interest to users of 'classic' MARK

Variance of a parameter average on the real scale

Postby lisark » Sun Jan 22, 2023 9:54 pm

Hi,

I am attempting to calculate the average and SE of four model estimates, subject to process variance alone. However, I am getting an error after the delta method. Could you please let me know (i) why the error occurs, and (ii) if my method is correct. At the end I have reproduced the calculation with total variance estimates without any problem.

My model has 16 Psi estimates, so I use Variance Components (random effects, intercept only) on all 16 Betas to obtain the 16 S-tilde which are subject only to process variance.

In order to calculate the average of the first 4 real Psi (subject only to process variance), I see that they can be calculated from the five logit-scaled S-tilde betas (labelled x1 to x5);
1st real psi = exp(x1+x2)/(1+exp(x1+x2))
2nd real psi = exp(x1+x3)/(1+exp(x1+x3))
3rd real psi = exp(x1+x4)/(1+exp(x1+x4))
4th real psi = exp(x1+x5)/(1+exp(x1+x5))

Then their average is;
Code: Select all
x1 = 1.480638; x2 = -0.897507; x3 = -1.033954; x4 = -0.673738; x5 = -0.582794
Average = (exp(x1+x2)/(1+exp(x1+x2)) + exp(x1+x3)/(1+exp(x1+x3)) + exp(x1+x4)/(1+exp(x1+x4)) + exp(x1+x5)/(1+exp(x1+x5)))/4 #0.6633982


And SE from the delta method is;
Code: Select all
vcov <- c(0.00208,    -0.50326,   -0.58330,    -0.60619,    -0.54785,
          -0.00221,     0.00922,     0.30596,     0.31018,     0.27697,
          -0.00218,     0.00241,     0.00673,    0.35633,     0.32541,
          -0.00229,     0.00247,     0.00242,     0.00686,    0.32754,
          -0.00211,     0.00224,     0.00225,     0.00229,     0.00712)

sigma <-matrix(vcov,nrow=5,ncol=5,byrow=TRUE)

seAv=deltamethod(~(exp(x1+x2)/(1+exp(x1+x2)) + exp(x1+x3)/(1+exp(x1+x3)) + exp(x1+x4)/(1+exp(x1+x4)) + exp(x1+x5)/(1+exp(x1+x5)))/4, c(x1,x2,x3,x4,x5), sigma, ses=T)

#Warning message: In sqrt(diag(new.covar)) : NaNs produced

CHECK WITH TOTAL VARIANCE ESTIMATES
MLE beta estimates:
Code: Select all
x1=1.490899; x2=-0.9146954; x3= -1.049563; x4=-0.6876908; x5=-0.5955833


Code: Select all
Av = (exp(x1+x2)/(1+exp(x1+x2)) + exp(x1+x3)/(1+exp(x1+x3)) + exp(x1+x4)/(1+exp(x1+x4)) + exp(x1+x5)/(1+exp(x1+x5)))/4 # 0.6623542


Code: Select all
vcov <- c(0.002135057, -0.002263782, -0.002240721, -0.002351096, -0.002162478,
-0.002263782,  0.009404006,  0.002471348,  0.002529644, 0.002297931,
-0.002240721,  0.002471348, 0.006842304,  0.002483221,  0.002307514,
-0.002351096,  0.002529644,  0.002483221,  0.006971759, 0.002345244,
-0.002162478,  0.002297931,  0.002307514,  0.002345244,  0.007232530)

sigma <-matrix(vcov,nrow=5,ncol=5,byrow=TRUE)

seAv=deltamethod(~(exp(x1+x2)/(1+exp(x1+x2)) + exp(x1+x3)/(1+exp(x1+x3)) + exp(x1+x4)/(1+exp(x1+x4)) + exp(x1+x5)/(1+exp(x1+x5)))/4, c(x1,x2,x3,x4,x5),sigma,ses=T) # 0.008148399
lisark
 
Posts: 4
Joined: Thu Jan 19, 2023 9:33 pm

Re: Variance of a parameter average on the real scale

Postby jlaake » Thu Feb 02, 2023 12:00 pm

I think you want to accomplish this by specifying a design argument to var.components or var.components.reml to take averages. I didn't look to see what your problem was to cause error.
jlaake
 
Posts: 1417
Joined: Fri May 12, 2006 12:50 pm
Location: Escondido, CA

Re: Variance of a parameter average on the real scale

Postby lisark » Mon Feb 13, 2023 9:53 pm

Thank you, var.components and var.components.reml gave much more sensible results.

My motivation for variance components with the whole dataset, then the delta method for a subset of the S-tilde was because of the advice to use at least 10-15 estimates for variance components analysis. The RMark dipper example uses only 6 and I have only 3 to 5, so is there a reason why var.components and var.components.reml are still appropriate?

I still think the delta method on the S-tilde should have worked ok, though obviously it didn't. Is there a reason why this approach should not be used?

Lastly, when var.components and var.components.reml yield similar but different results, should I prefer that from var.components.reml?
lisark
 
Posts: 4
Joined: Thu Jan 19, 2023 9:33 pm

Re: Variance of a parameter average on the real scale

Postby jlaake » Mon Feb 13, 2023 10:10 pm

Not sure i follow your logic. The stilde estimates were derived from var.components weren't they? So using those with delta method doesn't avoid concern about number of times being used to estimate process variance. My opinion that real is better but others may have thoughts.

Jeff
jlaake
 
Posts: 1417
Joined: Fri May 12, 2006 12:50 pm
Location: Escondido, CA

Re: Variance of a parameter average on the real scale

Postby lisark » Mon Feb 20, 2023 10:58 pm

Ok thank you. I guess having so few estimates is just a limitation in the data.

I was also thinking about finding a robust estimate of the median, rather than the mean.
Is it possible to do this with var.components or var.components.reml?
lisark
 
Posts: 4
Joined: Thu Jan 19, 2023 9:33 pm

Re: Variance of a parameter average on the real scale

Postby jlaake » Tue Feb 21, 2023 12:19 am

No there is not.
jlaake
 
Posts: 1417
Joined: Fri May 12, 2006 12:50 pm
Location: Escondido, CA


Return to RMark

Who is online

Users browsing this forum: Google [Bot] and 8 guests