question about statistical functions

[Mischkowski, Dominik (NIH/NCCIH) [F]]

Hi,

 

My name is Dominik Mischkowski, and I am postdoc with Lauren Atlas. We are currently working on a paper, indirectly comparing results using the functions glmfit_multilevel.m and mediation.m from the CAN lab statistics toolboxes.

 

We have noticed that we are getting diverging results for an IV predicting a DV when using the glmfit_multilevel and mediation functions, using the weighting and bootstrapping options. Here are the functions we are using:

 

glmfit_multilevel(DV, IV,[], 'names', {'Int', 'IV'}, 'weighted','boot','nresample',5000);

mediation(IV, DV, MEDIATOR, 'plots', 'verbose', 'boot','bootsamples',5000);

 

I also attached our scripts of these functions, so you know what exactly we are working with.

 

When running these functions, I find that the beta estimates of the IV-DV effects (in the mediation function, the c-path) are the same, but the estimates for the standard error of the beta (and consequently, the t-test) are different. Any idea what may the reason? Any thoughts on how to adjusting the functions so they give the same output when running a bootstrapping test to estimate MLM level coefficients?

 

Thank you very much!

-Dom

 

 

 

 

[Wani (Choong-Wan) Woo]

Hi Dom, 

Are they different a lot? or just a little bit different? Given that it's using bootstrap tests, which has a random component, it should be different every time in theory. 

Thanks,

Wani

--
You received this message because you are subscribed to the Google Groups "WagerlabTools" group.
To unsubscribe from this group and stop receiving emails from it, send an email to wagerlabtools+unsubscribe@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

--

Choong-Wan (Wani) Woo, PhD
Lab | Code | Google Scholar

Lab director/ Assistant professor

Computational Cognitive Affective Neuroscience Laboratory (Cocoan lab)

IBS Center for Neuroscience Imaging Research

Department of Biomedical Engineering

Sungkyunkwan University

Email: wan...@skku.edu 

          (or choong...@gmail.com)

Phone (office): +82 31-229-4363

[Mischkowski, Dominik (NIH/NCCIH) [F]]

Hi Wani,

 

Thanks for your quick response. No, it does not look like the differences are coming from an element of randomness introduced by the bootstrapping test.

 

Yes, they are substantially different. Not super, but substantially. We are not talking about slight differences in decimals here. If you want I can send you some output as well.

 

-Dom

 

 

 

From: "Wani (Choong-Wan) Woo" <choong...@gmail.com>
Date: Monday, May 15, 2017 at 8:48 PM
To: "Mischkowski, Dominik (NIH/NCCIH) [F]" <dominik.m...@nih.gov>
Cc: "wagerl...@googlegroups.com" <wagerl...@googlegroups.com>
Subject: Re: [WagerlabTools] question about statistical functions

 

Hi Dom, 

 

Are they different a lot? or just a little bit different? Given that it's using bootstrap tests, which has a random component, it should be different every time in theory. 

 

Thanks,

Wani

On Tue, May 16, 2017 at 4:17 AM, Mischkowski, Dominik (NIH/NCCIH) [F] <dominik.m...@nih.gov> wrote:

Hi,

 

My name is Dominik Mischkowski, and I am postdoc with Lauren Atlas. We are currently working on a paper, indirectly comparing results using the functions glmfit_multilevel.m and mediation.m from the CAN lab statistics toolboxes.

 

We have noticed that we are getting diverging results for an IV predicting a DV when using the glmfit_multilevel and mediation functions, using the weighting and bootstrapping options. Here are the functions we are using:

 

glmfit_multilevel(DV, IV,[], 'names', {'Int', 'IV'}, 'weighted','boot','nresample',5000);

mediation(IV, DV, MEDIATOR, 'plots', 'verbose', 'boot','bootsamples',5000);

 

I also attached our scripts of these functions, so you know what exactly we are working with.

 

When running these functions, I find that the beta estimates of the IV-DV effects (in the mediation function, the c-path) are the same, but the estimates for the standard error of the beta (and consequently, the t-test) are different. Any idea what may the reason? Any thoughts on how to adjusting the functions so they give the same output when running a bootstrapping test to estimate MLM level coefficients?

 

Thank you very much!

-Dom

 

 

 

 

--

You received this message because you are subscribed to the Google Groups "WagerlabTools" group.

To unsubscribe from this group and stop receiving emails from it, send an email to wagerlabtool...@googlegroups.com.

For more options, visit https://groups.google.com/d/optout.

[Wani (Choong-Wan) Woo]

Yep, it would be helpful if you can upload your variables and outputs here. So that people can run the exact code you ran. Thanks!!

- Wani

--

To unsubscribe from this group and stop receiving emails from it, send an email to wagerlabtools+unsubscribe@googlegroups.com.

For more options, visit https://groups.google.com/d/optout.

 

--

Choong-Wan (Wani) Woo, PhD
Lab | Code | Google Scholar

Lab director/ Assistant professor

Computational Cognitive Affective Neuroscience Laboratory (Cocoan lab)

IBS Center for Neuroscience Imaging Research

Department of Biomedical Engineering

Sungkyunkwan University

 

Email: wan...@skku.edu 

          (or choong...@gmail.com)

Phone (office): +82 31-229-4363

 

[Stephan Geuter]

Hi Dom,

I tried it on one dataset I had at hand. The difference was larger with 5000 bootstraps, but with 30k, it’s relatively small:

>> g.ste

ans =

   17.8048    0.3727

>> m.ste

ans =

    0.9850    0.0077    0.3597    0.3731    0.0579

5k:

>> g.ste

ans =

   17.9351    0.3749

>> m.ste

ans =

    0.9835    0.0078    0.3532    0.3697    0.0584

Best, Stephan

To unsubscribe from this group and stop receiving emails from it, send an email to wagerlabtool...@googlegroups.com.

[Mischkowski, Dominik (NIH/NCCIH) [F]]

Hi Stephan and Wani,

 

Thank you, that was very helpful. I reran bootstrapping analyses with more samples, and I got to the same conclusion. So far, looks like I underestimated the precision of a 5000 samples test.

 

But even with a smaller standard error, it looks like the t-test is less significant. Could this be due to a correction for more tests in the mediation function?

 

That is the output I got:

 

Mediation:

a b c' c ab

Coeff 0.80    259.68  249.98  484.87  190.63  

STE   0.04    70.04   71.31   38.84   46.81   

t (~N)    20.85   3.71    3.51    12.48   4.07    

Z 3.52    4.06    3.44    3.64    3.80    

p 0.0004  0.0000  0.0006  0.0003  0.0001

 

glmfit_multilevel

 

Int IV  

Coeff 6206.04 484.87  

STE   226.38  37.46   

t 30.09   13.56   

Z 3.49    3.47    

p 0.0002  0.0003  

 

But when I ran 30k, the standard error is almost identical, but the t-value is still off:

 

Mediation:

                a              b              c'             c              ab          

Coeff      0.80        259.68   249.98   484.87   190.63  

STE         0.04        70.02     70.69     38.19     46.85    

t (~N)     20.97     3.71        3.54        12.69     4.07       

Z              3.87        4.57        3.23        4.05        4.27       

p              0.0001   0.0000   0.0012   0.0001   0.0000  

 

glmfit_multilevel

                Int           IV           

Coeff      6206.04 484.87  

STE         227.56   38.23    

t              30.09     13.56    

Z              3.92        3.92       

p              0.0000   0.0000

 

 

For the sake of completion, I attached my analyses as a script and the accompanying variables.

 

I really appreciate you guys taking the time to respond.

 

Cheers,

-Dom

 

 

 

 

 

From: Stephan Geuter <Stephan...@Colorado.EDU>
Date: Monday, May 15, 2017 at 10:14 PM
To: "Mischkowski, Dominik (NIH/NCCIH) [F]" <dominik.m...@nih.gov>
Cc: "wagerl...@googlegroups.com" <wagerl...@googlegroups.com>, "Wani (Choong-Wan) Woo" <choong...@gmail.com>
Subject: Re: [WagerlabTools] question about statistical functions

 

Hi Dom,

I tried it on one dataset I had at hand. The difference was larger with 5000 bootstraps, but with 30k, it’s relatively small:

--

[Tor Wager]

Hi all,

No, it should not be due to any correction.  If you want to evaluate it, try it with the unweighted option first, and verify that you get the same results from the GLM. And I assume it’s path c you’re talking about.

You can evaluate bias and variance with some simulations, re-setting the RNG seed every time. There are some tools in the mediation toolbox for generating simulated data.  In any one dataset results will depend on the particular sample you have. At least 10K random datasets are needed to evaluate general properties.

Cheers,

Tor

<analyses.m><variables.mat>

[Mischkowski, Dominik (NIH/NCCIH) [F]]

Hi Tor et al.,

 

Yeah, right after I sent out the comment about multiple comparisons, I realized it was dumb. As far as I understand it (though maybe I am wrong), correcting for multiple comparisons should affect the dfs of t-value (and thus the p-value), not the t-value itself.

 

In any case, thank you, that was very helpful, too. I am going to run some simulations with real and simulated data.

 

Everybody have a great night! I appreciate y’all taking the time to respond.

 

-Dom

[Mischkowski, Dominik (NIH/NCCIH) [F]]

Hi everybody,

 

This is most likely my final email, closing this email thread.

 

Just to give some peace of mind to everybody involved – if needed: it seems that the weighting function makes the difference for calculating the t-statistic. When I run glmfit_multilevel and mediation using bootstrapping without weighting, the estimates for betas, their standard errors, and the t-values are essentially the same.

 

Thanks to Tor for the hint!

 

Everybody enjoy the rest of the week.

 

-Dom

[Mischkowski, Dominik (NIH/NCCIH) [F]]

Hi everybody,

 

Never say never. I apologize because I need to reopen this email thread. Unfortunately, an additional problem arouse, related to the issue raised previously.

 

I’d really appreciate anybody’s input of any kind here. I tried to figure out the problem myself – by debugging my functions – but my (still) limited capabilities in MatLab and my lack of understanding of the deep structure of these functions prohibited me from doing so.

 

Here is the background: When I ran a couple of analyses comparing the weighted with the non-weighted option in the glmfit_multilevel.m function (while also using 30K bootstrapping), I came across analysis output that did not make sense to me. I promise not to take up too much of your time here, but I think the issue may be of interest for other people using this function, too.

 

Before I get started: I attached (a) an analysis script depicting functions with their parameters (also includes complete output), (b) the variables used in my analyses, and (c) both the glmfit_multilevel and _general functions that I am using. I downloaded these functions about 1.5 years ago from the CAN lab website (thanks!) and I did not temper with them, but just in case.

 

So, here is the problem:

 

When I ran MLM analyses with weighting, I got output that I did not understand:

 

[stats_weighted = glmfit_multilevel(X9, X2,[], 'names', {'Int', 'IV'},'weighted','boot','nresample',30000);]

 

2nd-level B01

    Int IV  

Coeff   6819.00 13.89   

STE 502.32  8.77    

t   13.42   0.84    

Z   4.08    1.37    

p   0.0000  0.0851

 

In this output, the t-value (bold) of the effect of the IV did not make sense to me, because of two reasons: 1) the coefficient > its STE, which is a good indicator that the t-value should have been > 1 – after all, t is basically calculated by setting a coefficient (e.g., beta) in relation to its STE, and 2) the z-value was bigger (!) than t, which theoretically does not make sense. At any given p-value, t needs to approximate z from the upper bound, since t models uncertainty from sampling error in addition to z.

 

When I ran the same analyses without weighting, t-values, z-values and p-values looked OK to me (see output in attached analysis_script).

 

As a side question: Is it possible to run these analyses two-sided instead of one-sided? At least that is how I understand the z-value – its relation to the p-value only makes sense when z is tested one-sided. I understand that there are philosophical differences in whether to run one- vs. two-sided tests, but at least in my analyses, I can’t assume that effects will be only in one direction.

 

I hope with all these questions I did not miss the obvious. Sometimes, putting the right parameters into the function does wonders, after all!

 

Anyways, I will be eternally thankful for any input and will pay in alcohol at the next conference.

 

Everybody have a good night!

-Dom

[Mischkowski, Dominik (NIH/NCCIH) [F]]

Hi,

 

I apologize: Maybe I was too complicated in my previous email, and thus did not receive a response. So, I am trying to do better:

 

Lauren and I right now are trying to parameterize glmfit_multilevel’ so it provides the same analyses as the mediation function with parameters: bootstrapping and weighting. But we are having trouble with the output of the glmfit_multilevel function:

 

The t-test is off (<1, but still marginal significant):

 

(stats = glmfit_multilevel(X9, X2,[], 'names', {'Int', 'pain'}, 'analysisname', 'SCR AUC (0-4s) by pain','weighted','boot','nresample',30000);)

                Int           X             

Coeff      6819.00 13.89    

STE         505.39   8.69       

t              13.42     0.84       

Z              4.08        1.41       

p              0.0000   0.0791  

 

In contrast: When we run mediation.m using bootstrapping and weighting (should be equivalent with the previous test), we get similar coefficients and STEs (see path c), but the t-test looks much more sensible:

 

(I added an artificial mediator variable to be able to calculate the overall c-path – disregard the remaining paths)

([paths, stats] = mediation(X2, X9, X1, 'plots', 'verbose', 'names', {'X', 'Y', 'Mediator'}, 'plots','boot','bootsamples',30000);)

Multi-level model

                a              b              c'             c              ab          

Coeff      0.92        -60.26    72.59     13.89     -48.31   

STE         0.05        21.71     20.08     8.69        17.26    

t (~N)     17.72     -2.78      3.62        1.60        -2.80     

Z              4.20        -2.70      4.19        1.78        -2.76     

p              0.0000   0.0069   0.0000   0.0743   0.0058  

 

 

For comparison: When we run glmfit_multilevel, using bootstrapping but no weighting, we also get sensible results:

(stats = glmfit_multilevel(X9, X2,[], 'names', {'Int', 'X'},'boot','nresample',30000);)

 

2nd-level B01

                Int           X             

Coeff      6827.00 37.32    

STE         506.95   19.96    

t              13.43     1.85       

Z              4.09        1.83       

p              0.0000   0.0339

 

So, the problem really is with the weighting parameter in the glmfit_multilevel function, but I have trouble figuring out why. I hope this makes sense. Any advice is appreciated!

 

Cheers,

-Dominik