I'm trying to calculate an effect size for MIXED according to Selya et al 2012 paper. For this
method, it's important to obtain the residual variance of the given model. Now for a simple model without a repeated covariance structure, the residual variance of the model is reported in the output, like below:
However, once there is a covariance structure added, like ARMA11, I am not quite sure how we're supposed to obtain the residual variance for the model. Here is a sample output: The same goes for the autoregressive model: Does anyone know how are the residual variances obtained for these more complex models? Or perhaps, alternatively, if there is a different method of calculating an effect size. The issue is the same when attempting a calculation of the ICC. E
CONFIDENTIALITY NOTICE: This e-mail communication and any attachments may contain confidential and privileged information for the use of the designated recipients named above. If you are not the intended recipient, you are hereby notified that you have received this communication in error and that any review, disclosure, dissemination, distribution or copying of it or its contents is prohibited. If you have received this communication in error, please notify me immediately by replying to this message and destroy all copies of this communication and any attachments. Thank you.
=====================
To manage your subscription to SPSSX-L, send a message to
[hidden email] (not to SPSSX-L), with no body text except the
command. To leave the list, send the command
SIGNOFF SPSSX-L
For a list of commands to manage subscriptions, send the command
INFO REFCARD
|
Emil – You will probably have better chance of getting this answered on the Multilevel List (https://www.jiscmail.ac.uk/cgi-bin/webadmin?A0=MULTILEVEL). Dr. Hedeker responds to that list. Maybe contact the author(s) directly, too? peter From: SPSSX(r) Discussion [mailto:[hidden email]] On Behalf Of Rudobeck, Emil (LLU) I'm trying to calculate an effect size for MIXED according to Selya et al 2012 paper. For this method, it's important to obtain the residual variance of the given model. Now for a simple model without a repeated covariance structure, the residual variance of the model is reported in the output, like below: CONFIDENTIALITY NOTICE: This e-mail communication and any attachments may contain confidential and privileged information for the use of the designated recipients named above. If you are not the intended recipient, you are hereby notified that you have received this communication in error and that any review, disclosure, dissemination, distribution or copying of it or its contents is prohibited. If you have received this communication in error, please notify me immediately by replying to this message and destroy all copies of this communication and any attachments. Thank you. ===================== To manage your subscription to SPSSX-L, send a message to [hidden email] (not to SPSSX-L), with no body text except the command. To leave the list, send the command SIGNOFF SPSSX-L For a list of commands to manage subscriptions, send the command INFO REFCARD |
I can't see the output from SPSS-L, but the residual variance should be provided when specifying ARMA11 or AR1 on the REPEATED statement. Ryan On Wed, Jul 13, 2016 at 1:03 PM, Peter Link <[hidden email]> wrote:
|
Here is an example I found in the SPSS-L archives a long time ago where I generated data that arose from a first-order autoregressive residual covariance structure by taking advantage of the cholesky decomposition and then fit the model using MIXED. If you run the code you will see the estimated residual covariance is provided. I don't have time at the moment to discuss effect sizes but I figured the example below might help you get started. Unless someone else responds, I will try to address your effect size question at a later date. Ryan -- *Generate Data for Mixed Model with AR1 specification. set seed 65923454. new file. inp pro. compute subject=-99. compute time = -99. compute x1 = -99. compute x2 = -99. compute x3 = -99. compute e1 = -99. compute e2 = -99. compute e3 = -99. compute sigma = 1. compute rho = 0.50. compute a11 = 1. compute a21 = rho. compute a31 = rho**2. compute a22 = sqrt(1 - rho**2). compute a32 = rho*sqrt(1 - rho**2). compute a33 = sqrt(1 - rho**2). leave subject to a33. loop subject= 1 to 1000. compute x1 = rv.normal(0,1). compute x2 = rv.normal(0,1). compute x3 = rv.normal(0,1). compute e1 = sigma * a11*x1. compute e2 = sigma * (a21*x1 + a22*x2). compute e3 = sigma * (a31*x1 + a32*x2 + a33*x3). loop time = 1 to 3. compute y = 1.5 + 1.2*(time=1) + 0.9*(time=2) + e1*(time=1) + e2*(time=2) + e3*(time=3). end case. end loop. end loop. end file. end inp pro. exe. delete variables x1 x2 x3 sigma rho a11 a21 a31 a22 a32 a33 e1 e2 e3. MIXED y BY time /FIXED=time | SSTYPE(3) /METHOD=REML /PRINT=R SOLUTION /REPEATED=time | SUBJECT(subject) COVTYPE(AR1). On Wed, Jul 13, 2016 at 1:28 PM, Ryan Black <[hidden email]> wrote:
|
OP and others: Although time is against me, I figured I'd estimate the pseudo R-squared estimate associated with the fixed effect of time from the model with an AR1 specification I posted earlier today. For this model, pseudo R-squared can be calculated as: (unconditional residual variance - conditional residual variance) / unconditional residual variance OR 1 - (conditional residual variance / unconditional residual variance) Below is the code. Feel free to write back if you have questions. I'll try to respond as time permits. HTH. Ryan -- *Generate Data for Mixed Model with AR1 specification. set seed 65923454. new file. inp pro. compute subject=-99. compute time = -99. compute x1 = -99. compute x2 = -99. compute x3 = -99. compute e1 = -99. compute e2 = -99. compute e3 = -99. compute sigma = 1. compute rho = 0.50. compute a11 = 1. compute a21 = rho. compute a31 = rho**2. compute a22 = sqrt(1 - rho**2). compute a32 = rho*sqrt(1 - rho**2). compute a33 = sqrt(1 - rho**2). leave subject to a33. loop subject= 1 to 1000. compute x1 = rv.normal(0,1). compute x2 = rv.normal(0,1). compute x3 = rv.normal(0,1). compute e1 = sigma * a11*x1. compute e2 = sigma * (a21*x1 + a22*x2). compute e3 = sigma * (a31*x1 + a32*x2 + a33*x3). loop time = 1 to 3. compute y = 1.5 + 1.2*(time=1) + 0.9*(time=2) + e1*(time=1) + e2*(time=2) + e3*(time=3). end case. end loop. end loop. end file. end inp pro. exe. delete variables x1 x2 x3 sigma rho a11 a21 a31 a22 a32 a33 e1 e2 e3. MIXED y BY time /FIXED= | SSTYPE(3) /METHOD=REML /PRINT=R SOLUTION /REPEATED=time | SUBJECT(subject) COVTYPE(AR1). MIXED y BY time /FIXED=time | SSTYPE(3) /METHOD=REML /PRINT=R SOLUTION /REPEATED=time | SUBJECT(subject) COVTYPE(AR1). compute pseudo_R_squared_1 = (1.319641 - 1.035615) / 1.319641. compute pseudo_R_squared_2 = 1 - (1.035615 / 1.319641). execute. On Wed, Jul 13, 2016 at 1:49 PM, Ryan Black <[hidden email]> wrote:
|
Emil, Please read this post carefully. 1. I'm generating data which conform to the specified MIXED model. Data simulation on a forum helps facilitate communication among all of us without having to look at a particular dataset. I can explain how the data are generated in more detail but given your previous comment I suspect it will not help you get any closer to what you need to know. For now, please just RUN the code that generates the data and LOOK at the dataset. It should resemble something similar to what you are observing in your dataset. (I don't have time to review your dataset.) 2. ARMA11 and AR1 both assume homogeneous variances as stated on the IBM website: Specifically, the website states the following: "ARMA(1,1). This is a first-order autoregressive moving average structure. It has homogenous variances. The correlation between two elements is equal to f*r for adjacent elements, f*(r2) for elements separated by a third, and so on. r and f are the autoregressive and moving average parameters, respectively, and their values are constrained to lie between –1 and 1, inclusive." In other words, this structure produces an estimate of the error variance (aka residual variance) as well as estimates of autoregressive and moving parameters. Therefore, my illustration is applicable to ARMA11. Change COVTYPE(AR1) to COVTYPE(ARMA11) in the MIXED code I provided previously, obtain the unconditional residual variance (without time effect) and the conditional residual variance (with the time effect) and plug those estimates into one of the two formulas I provided previously to obtain an estimate of the pseudo R-squared associated with inclusion of the fixed effect of time. 3. Now, if you have changed your mind in that you no longer want to specify ARMA11 as the residual covariance type, then please update SPSS-L with the specific residual covariance type you want to specify. Ryan On Thu, Jul 14, 2016 at 2:56 AM, Rudobeck, Emil (LLU) <[hidden email]> wrote:
|
Administrator
|
In reply to this post by Ryan
Some enhancements to that really ugly non general simulation code ;-)
-- DEFINE SimulateAR1_RBlack (Sigma !TOKENS(1)/Rho !TOKENS(1) /N !TOKENS(1) /P !TOKENS(1) /Const !TOKENS(1) /Beta !CMDEND ). /* Create Canvas */. MATRIX. SAVE MAKE(!N,!P,0)/OUTFILE * /VARIABLES X1 TO !CONCAT(X,!P). END MATRIX. /* Populate with Normal(0,1) and ID */. DO REPEAT X=X1 TO !CONCAT(X,!P). COMPUTE X=RV.NORMAL(0,1). END REPEAT. COMPUTE ID=$CASENUM. MATRIX. GET X /FILE * /VARIABLES X1 TO !CONCAT(X,!P). GET ID /FILE * /VARIABLES ID . /*Build AR covariance structure */. COMPUTE r=MAKE(!p,!p,!Rho). CALL SETDIAG(r,1). LOOP #=3 TO !P. + LOOP ##=1 TO #-1. + COMPUTE r(#,##)=r(#,##)*r(#-1,##). + COMPUTE r(##,#)=r(#,##). + END LOOP. END LOOP. /* PRINT r . /* PRINT CHOL(R). /* Build time variable */. COMPUTE Time=MAKE(!P,1,1). LOOP #=1 TO !P. + COMPUTE Time(#,1)=#. END LOOP. SAVE ({KRONEKER(MAKE(!N,1,1),Time), KRONEKER({ID, X * CHOL(r) * !Sigma},MAKE(!P,1,1)) }) /OUTFILE * /VARIABLES Time ID E1 TO !CONCAT(E,!P). END MATRIX. VECTOR E = E1 TO !CONCAT(E,!P). COMPUTE Y=SUM(!Const, E(Time)). DO REPEAT T=1 TO !P / B=!Beta. COMPUTE Y=SUM(Y, (Time EQ T) * B). END REPEAT . !ENDDEFINE. SimulateAR1_RBlack Sigma=1 Rho=.5 N=100 P=3 Const=1.5 Beta= 1.2 .9 0 .
Please reply to the list and not to my personal email.
Those desiring my consulting or training services please feel free to email me. --- "Nolite dare sanctum canibus neque mittatis margaritas vestras ante porcos ne forte conculcent eas pedibus suis." Cum es damnatorum possederunt porcos iens ut salire off sanguinum cliff in abyssum?" |
Administrator
|
Don't hold back, David. Tell us what you really thought of Ryan's simulation code. ;-)
Regarding that macro call... SimulateAR1_RBlack Sigma=1 Rho=.5 N=100 P=3 Const=1.5 Beta= 1.2 .9 0 . ...note that the number of arguments after 'Beta =' is 3, the same value that was assigned to P (the number of predictor variables). HTH.
--
Bruce Weaver bweaver@lakeheadu.ca http://sites.google.com/a/lakeheadu.ca/bweaver/ "When all else fails, RTFM." PLEASE NOTE THE FOLLOWING: 1. My Hotmail account is not monitored regularly. To send me an e-mail, please use the address shown above. 2. The SPSSX Discussion forum on Nabble is no longer linked to the SPSSX-L listserv administered by UGA (https://listserv.uga.edu/). |
Administrator
|
As soon as you start doing Matrix algebra like that you know yer in trouble.
Try generalizing that to 4 time periods -Once upon a time [LONG AGO]- I sat down and worked out explicit formulas for each element of the inverse of a 3 x 3 Matrix. Then I got clever and tried it for a 4x4 and it took up a poster sized sheet of paper with rather small print ;-)
Please reply to the list and not to my personal email.
Those desiring my consulting or training services please feel free to email me. --- "Nolite dare sanctum canibus neque mittatis margaritas vestras ante porcos ne forte conculcent eas pedibus suis." Cum es damnatorum possederunt porcos iens ut salire off sanguinum cliff in abyssum?" |
In reply to this post by David Marso
SPSS-L: We can construct e ~ MVN(0, SIGMA) by computing e = A*X where A is the lower triangular matrix of a Cholesky decomposition of SIGMA and X ~ MVN(0, I(k)). If _ _ with _ _ then the following values of the matrix A can be shown to hold: a11 = 1
e1 = sigma * a11*x1; Hope this helps provide clarification to the OP, although I do not think it's necessary given the OP's request to obtain a standardized effect size. Thanks to David for providing more efficient and generalizable code. Ryan On Thu, Jul 14, 2016 at 10:51 AM, David Marso <[hidden email]> wrote: Some enhancements to that really ugly non general simulation code ;-) |
Small correction in the message I just posted. I wrote this too fast because I'm trying to get back to work. We can construct e ~ MVN(0, SIGMA) by computing e = A*X where A is the lower triangular matrix of a Cholesky decomposition of SIGMA and X ~ MVN(0, I(k)). If _ _ with _ _ then the following values of the matrix A can be shown to hold: a11 = 1
e1 = sigma * a11*x1 e2 = sigma * (a21*x1 + a22*x2) e3 = sigma * (a31*x1 + a32*x2 + a33*x3) From there, one can see how I constructed the "ugly" SPSS code. :-) Now it's time for me to get back to work! Ryan On Fri, Jul 15, 2016 at 8:19 AM, Ryan Black <[hidden email]> wrote:
|
In reply to this post by Ryan
I think my last message was rejected because of attachments. I don't know what types of attachments are allowed, if any, but here is a link to the files:
https://drive.google.com/folderview?id=0B0AZLzY6tv15azBIaDU1bjZyU28&usp=sharing
And here is the message again: Hi Ryan, Thanks for the details. I am not that great with SPSS syntax. I think in your code you're just generating a dataset from random numbers and calculating based on those values, but I am not sure about the details, like why do you set x, e values = -99, unless it's just random designation. Also, rho = 0.5 presumably would be substituted with the actual rho value from unconditional the model. This example would be easier to understand if we get rid of the random number syntax. To that end, I have attached a simplified dataset along with the MIXED syntax for both AR1 and ARMA11. It would be helpful to see how your syntax for effect size would integrate with such a sample dataset to see exactly what is being calculated. Some questions: 1. A lot of times I encounter heterogeneous structure as the best fit. I'm assuming this effect size method applies not only to AR1, but to ARH1 as well. 2. My original issue was and still is with ARMA11. Are you aware of the methods for calculating the residual variance and/or the effect size for ARMA11 structures? In addition to rho, there is also phi for ARMA11, so AR1 effect size syntax probably can't be used for it. Thanks, Emil From: SPSSX(r) Discussion [[hidden email]] on behalf of Ryan Black [[hidden email]]
Sent: Wednesday, July 13, 2016 5:54 PM To: [hidden email] Subject: Re: [BULK] [SPSSX-L] Mixed Model Effect Size OP and others:
Although time is against me, I figured I'd estimate the pseudo R-squared estimate associated with the fixed effect of time from the model with an AR1 specification I posted earlier today.
For this model, pseudo R-squared can be calculated as:
(unconditional residual variance - conditional residual variance) / unconditional residual variance
OR
1 - (conditional residual variance / unconditional residual variance)
Below is the code.
Feel free to write back if you have questions. I'll try to respond as time permits.
HTH.
Ryan
--
*Generate Data for Mixed Model with AR1 specification.
set seed 65923454. new file.
inp pro. compute subject=-99.
compute time = -99. compute x1 = -99. compute x2 = -99. compute x3 = -99. compute e1 = -99. compute e2 = -99. compute e3 = -99. compute sigma = 1. compute rho = 0.50. compute a11 = 1. compute a21 = rho. compute a31 = rho**2. compute a22 = sqrt(1 - rho**2). compute a32 = rho*sqrt(1 - rho**2). compute a33 = sqrt(1 - rho**2). leave subject to a33.
loop subject= 1 to 1000.
compute x1 = rv.normal(0,1).
compute x2 = rv.normal(0,1). compute x3 = rv.normal(0,1). compute e1 = sigma * a11*x1.
compute e2 = sigma * (a21*x1 + a22*x2). compute e3 = sigma * (a31*x1 + a32*x2 + a33*x3). loop time = 1 to 3.
compute y = 1.5 + 1.2*(time=1) + 0.9*(time=2) + e1*(time=1) +
e2*(time=2) + e3*(time=3). end case. end loop. end loop. end file. end inp pro. exe. delete variables x1 x2 x3 sigma rho a11 a21 a31 a22 a32 a33 e1 e2 e3.
MIXED y BY time
/FIXED= | SSTYPE(3) /METHOD=REML /PRINT=R SOLUTION /REPEATED=time | SUBJECT(subject) COVTYPE(AR1). MIXED y BY time
/FIXED=time | SSTYPE(3) /METHOD=REML /PRINT=R SOLUTION /REPEATED=time | SUBJECT(subject) COVTYPE(AR1). compute pseudo_R_squared_1 = (1.319641 - 1.035615) / 1.319641.
compute pseudo_R_squared_2 = 1 - (1.035615 / 1.319641). execute. On Wed, Jul 13, 2016 at 1:49 PM, Ryan Black
<ryan.andrew.black@...> wrote:
WARNING: Please be vigilant when opening emails that appear to be the least bit out of the ordinary, e.g. someone you usually don’t hear from, or attachments you usually don’t receive or didn’t expect, requests to click links or log into systems, etc. If you receive suspicious emails, please do not open attachments or links and immediately forward the suspicious email to [hidden email] and then delete the suspicious email. |
Never mind my last message. I just stumbled upon this very same thread in the SPSSX online forum and saw responses to my message that I never received via email (I'm using only
email for this communication and hadn't checked out the forum much)...
=====================
To manage your subscription to SPSSX-L, send a message to
[hidden email] (not to SPSSX-L), with no body text except the
command. To leave the list, send the command
SIGNOFF SPSSX-L
For a list of commands to manage subscriptions, send the command
INFO REFCARD
Alright, so I reviewed the code more closely. I had mistakenly assumed it's a syntax for some complex calculation of an effect size (probably because I saw sigma and rho there), but all that syntax is just for generating a Gaussian dataset. Correct me if I'm wrong, but it just boils down to this: the residual variance of AR1 and ARMA11 structures is equivalent to the constant variance across the repeated measures (diagonal). Is that it? If so, it would be good to see some sort of a validating source on this because this method ignores other errors, such as the autocorrelation errors between the measurements (rho). Maybe that's not important, but I can't be sure without a reference. I'm also not clear why the syntax ends with two pseudo-r square calculations (equivalent output). You can only get one r square from this set of two models, but let me know if I'm missing something. I know that AR1 is not heterogeneous. In my original message, that was clipped, I was referring to ARH1. But we can get to that after clarifying the effect size for the homogeneous structures. The type of structure is dataset dependent and I have found a best fit with ARMA11, AR1, and even AD1 and TP structures. Emil From: Rudobeck, Emil (LLU)
Sent: Friday, July 15, 2016 8:30 PM To: Ryan Black; [hidden email] Subject: RE: [BULK] [SPSSX-L] Mixed Model Effect Size I think my last message was rejected because of attachments. I don't know what types of attachments are allowed, if any, but here is a link to the files:
https://drive.google.com/folderview?id=0B0AZLzY6tv15azBIaDU1bjZyU28&usp=sharing
And here is the message again: Hi Ryan, Thanks for the details. I am not that great with SPSS syntax. I think in your code you're just generating a dataset from random numbers and calculating based on those values, but I am not sure about the details, like why do you set x, e values = -99, unless it's just random designation. Also, rho = 0.5 presumably would be substituted with the actual rho value from unconditional the model. This example would be easier to understand if we get rid of the random number syntax. To that end, I have attached a simplified dataset along with the MIXED syntax for both AR1 and ARMA11. It would be helpful to see how your syntax for effect size would integrate with such a sample dataset to see exactly what is being calculated. Some questions: 1. A lot of times I encounter heterogeneous structure as the best fit. I'm assuming this effect size method applies not only to AR1, but to ARH1 as well. 2. My original issue was and still is with ARMA11. Are you aware of the methods for calculating the residual variance and/or the effect size for ARMA11 structures? In addition to rho, there is also phi for ARMA11, so AR1 effect size syntax probably can't be used for it. Thanks, Emil From: SPSSX(r) Discussion [[hidden email]] on behalf of Ryan Black [[hidden email]]
Sent: Wednesday, July 13, 2016 5:54 PM To: [hidden email] Subject: Re: [BULK] [SPSSX-L] Mixed Model Effect Size OP and others:
Although time is against me, I figured I'd estimate the pseudo R-squared estimate associated with the fixed effect of time from the model with an AR1 specification I posted earlier today.
For this model, pseudo R-squared can be calculated as:
(unconditional residual variance - conditional residual variance) / unconditional residual variance
OR
1 - (conditional residual variance / unconditional residual variance)
Below is the code.
Feel free to write back if you have questions. I'll try to respond as time permits.
HTH.
Ryan
--
*Generate Data for Mixed Model with AR1 specification.
set seed 65923454. new file.
inp pro. compute subject=-99.
compute time = -99. compute x1 = -99. compute x2 = -99. compute x3 = -99. compute e1 = -99. compute e2 = -99. compute e3 = -99. compute sigma = 1. compute rho = 0.50. compute a11 = 1. compute a21 = rho. compute a31 = rho**2. compute a22 = sqrt(1 - rho**2). compute a32 = rho*sqrt(1 - rho**2). compute a33 = sqrt(1 - rho**2). leave subject to a33.
loop subject= 1 to 1000.
compute x1 = rv.normal(0,1).
compute x2 = rv.normal(0,1). compute x3 = rv.normal(0,1). compute e1 = sigma * a11*x1.
compute e2 = sigma * (a21*x1 + a22*x2). compute e3 = sigma * (a31*x1 + a32*x2 + a33*x3). loop time = 1 to 3.
compute y = 1.5 + 1.2*(time=1) + 0.9*(time=2) + e1*(time=1) +
e2*(time=2) + e3*(time=3). end case. end loop. end loop. end file. end inp pro. exe. delete variables x1 x2 x3 sigma rho a11 a21 a31 a22 a32 a33 e1 e2 e3.
MIXED y BY time
/FIXED= | SSTYPE(3) /METHOD=REML /PRINT=R SOLUTION /REPEATED=time | SUBJECT(subject) COVTYPE(AR1). MIXED y BY time
/FIXED=time | SSTYPE(3) /METHOD=REML /PRINT=R SOLUTION /REPEATED=time | SUBJECT(subject) COVTYPE(AR1). compute pseudo_R_squared_1 = (1.319641 - 1.035615) / 1.319641.
compute pseudo_R_squared_2 = 1 - (1.035615 / 1.319641). execute. On Wed, Jul 13, 2016 at 1:49 PM, Ryan Black
<ryan.andrew.black@...> wrote:
WARNING: Please be vigilant when opening emails that appear to be the least bit out of the ordinary, e.g. someone you usually don’t hear from, or attachments you usually don’t receive or didn’t expect, requests to click links or log into systems, etc. If you receive suspicious emails, please do not open attachments or links and immediately forward the suspicious email to [hidden email] and then delete the suspicious email. |
I'm going to bump this. If no one else here knows that much about repeated LMM effect size, I will start contacting authors whose books I use for reference quite often. Just wanted
to be sure it's not common knowledge first, which given by the lack of responses seems not to be. I will pursue the multilevel forum as well, as was suggested by Peter.
Ryan, I'm not sure if you don't know more about the LMM effect size or if you're probably busy. In either case, I appreciate the effort to help. Sorry that I had missed some of your responses before. I welcome any practical references for repeated LMM effect size calculations. Note that I did not find Raudenbush or Snijders (as referenced in this technote http://www-01.ibm.com/support/docview.wss?uid=swg21474983) useful to that end. E From: Rudobeck, Emil (LLU)
Sent: Sunday, July 17, 2016 12:49 AM To: Ryan Black; [hidden email] Subject: RE: [BULK] [SPSSX-L] Mixed Model Effect Size Never mind my last message. I just stumbled upon this very same thread in the SPSSX online forum and saw responses to my message that I never received via email (I'm using only email
for this communication and hadn't checked out the forum much)...
Alright, so I reviewed the code more closely. I had mistakenly assumed it's a syntax for some complex calculation of an effect size (probably because I saw sigma and rho there), but all that syntax is just for generating a Gaussian dataset. Correct me if I'm wrong, but it just boils down to this: the residual variance of AR1 and ARMA11 structures is equivalent to the constant variance across the repeated measures (diagonal). Is that it? If so, it would be good to see some sort of a validating source on this because this method ignores other errors, such as the autocorrelation errors between the measurements (rho). Maybe that's not important, but I can't be sure without a reference. I'm also not clear why the syntax ends with two pseudo-r square calculations (equivalent output). You can only get one r square from this set of two models, but let me know if I'm missing something. I know that AR1 is not heterogeneous. In my original message, that was clipped, I was referring to ARH1. But we can get to that after clarifying the effect size for the homogeneous structures. The type of structure is dataset dependent and I have found a best fit with ARMA11, AR1, and even AD1 and TP structures. Emil From: Rudobeck, Emil (LLU)
Sent: Friday, July 15, 2016 8:30 PM To: Ryan Black; [hidden email] Subject: RE: [BULK] [SPSSX-L] Mixed Model Effect Size I think my last message was rejected because of attachments. I don't know what types of attachments are allowed, if any, but here is a link to the files:
https://drive.google.com/folderview?id=0B0AZLzY6tv15azBIaDU1bjZyU28&usp=sharing
And here is the message again: Hi Ryan, Thanks for the details. I am not that great with SPSS syntax. I think in your code you're just generating a dataset from random numbers and calculating based on those values, but I am not sure about the details, like why do you set x, e values = -99, unless it's just random designation. Also, rho = 0.5 presumably would be substituted with the actual rho value from unconditional the model. This example would be easier to understand if we get rid of the random number syntax. To that end, I have attached a simplified dataset along with the MIXED syntax for both AR1 and ARMA11. It would be helpful to see how your syntax for effect size would integrate with such a sample dataset to see exactly what is being calculated. Some questions: 1. A lot of times I encounter heterogeneous structure as the best fit. I'm assuming this effect size method applies not only to AR1, but to ARH1 as well. 2. My original issue was and still is with ARMA11. Are you aware of the methods for calculating the residual variance and/or the effect size for ARMA11 structures? In addition to rho, there is also phi for ARMA11, so AR1 effect size syntax probably can't be used for it. Thanks, Emil From: SPSSX(r) Discussion [[hidden email]] on behalf of Ryan Black [[hidden email]]
Sent: Wednesday, July 13, 2016 5:54 PM To: [hidden email] Subject: Re: [BULK] [SPSSX-L] Mixed Model Effect Size OP and others:
Although time is against me, I figured I'd estimate the pseudo R-squared estimate associated with the fixed effect of time from the model with an AR1 specification I posted earlier today.
For this model, pseudo R-squared can be calculated as:
(unconditional residual variance - conditional residual variance) / unconditional residual variance
OR
1 - (conditional residual variance / unconditional residual variance)
Below is the code.
Feel free to write back if you have questions. I'll try to respond as time permits.
HTH.
Ryan
--
*Generate Data for Mixed Model with AR1 specification.
set seed 65923454. new file.
inp pro. compute subject=-99.
compute time = -99. compute x1 = -99. compute x2 = -99. compute x3 = -99. compute e1 = -99. compute e2 = -99. compute e3 = -99. compute sigma = 1. compute rho = 0.50. compute a11 = 1. compute a21 = rho. compute a31 = rho**2. compute a22 = sqrt(1 - rho**2). compute a32 = rho*sqrt(1 - rho**2). compute a33 = sqrt(1 - rho**2). leave subject to a33.
loop subject= 1 to 1000.
compute x1 = rv.normal(0,1).
compute x2 = rv.normal(0,1). compute x3 = rv.normal(0,1). compute e1 = sigma * a11*x1.
compute e2 = sigma * (a21*x1 + a22*x2). compute e3 = sigma * (a31*x1 + a32*x2 + a33*x3). loop time = 1 to 3.
compute y = 1.5 + 1.2*(time=1) + 0.9*(time=2) + e1*(time=1) +
e2*(time=2) + e3*(time=3). end case. end loop. end loop. end file. end inp pro. exe. delete variables x1 x2 x3 sigma rho a11 a21 a31 a22 a32 a33 e1 e2 e3.
MIXED y BY time
/FIXED= | SSTYPE(3) /METHOD=REML /PRINT=R SOLUTION /REPEATED=time | SUBJECT(subject) COVTYPE(AR1). MIXED y BY time
/FIXED=time | SSTYPE(3) /METHOD=REML /PRINT=R SOLUTION /REPEATED=time | SUBJECT(subject) COVTYPE(AR1). compute pseudo_R_squared_1 = (1.319641 - 1.035615) / 1.319641.
compute pseudo_R_squared_2 = 1 - (1.035615 / 1.319641). execute. On Wed, Jul 13, 2016 at 1:49 PM, Ryan Black
<ryan.andrew.black@...> wrote:
WARNING: Please be vigilant when opening emails that appear to be the least bit out of the ordinary, e.g. someone you usually don’t hear from, or attachments you usually don’t receive or didn’t expect, requests to click links or log into systems, etc. If you receive suspicious emails, please do not open attachments or links and immediately forward the suspicious email to [hidden email] and then delete the suspicious email. |
I want to give this another try. After reading more and exploring SAS, I wanted to discuss some specific options.
=====================
To manage your subscription to SPSSX-L, send a message to
[hidden email] (not to SPSSX-L), with no body text except the
command. To leave the list, send the command
SIGNOFF SPSSX-L
For a list of commands to manage subscriptions, send the command
INFO REFCARD
1. Does anyone have experience with converting LMM SEs into SD (SD=SE*sqrt[df+1]) and reporting Cohen's d, at least for the contrasts? Any specific caveats? 2. Another option I've been contemplating is to extract predicted values (PRED) from LMM and calculate whatever effect size I need. As an example, I could run a t-test on PRED values for a given contrast of interest and then report Cohen's d. This could be also extended to running an ANOVA on all predicted value categories and reporting partial eta squared as well. The above solutions aren't ideal, but I have not been able to find any other options. I think Ryan's suggestion of pseudo-R squared is based on the methods proposed by Selya et al that I referred to before. Those methods are not designed for complex covariance structures like AR1, which doesn't just have the autocorrelation coefficient (rho), but also the diagonal residual (2 sources of errors). ARMA11 has 3 sources of errors, from what I understand. I also got in touch with Selya and verified that the methods there don't apply to the complex covariance structures. I did post in the MLM forum and did not get any responses. For whatever reason, LMM effect size appears to have been overlooked by statisticians, especially given that mixed models have been around for over 150 years now, apparently first formulated by astronomer George Airy. We need someone like Cohen to come up with a standardized effect size for mixed models. From: Rudobeck, Emil (LLU)
Sent: Thursday, July 21, 2016 8:59 PM To: [hidden email] Subject: RE: [BULK] [SPSSX-L] Mixed Model Effect Size I'm going to bump this. If no one else here knows that much about repeated LMM effect size, I will start contacting authors whose books I use for reference quite often. Just wanted
to be sure it's not common knowledge first, which given by the lack of responses seems not to be. I will pursue the multilevel forum as well, as was suggested by Peter.
Ryan, I'm not sure if you don't know more about the LMM effect size or if you're probably busy. In either case, I appreciate the effort to help. Sorry that I had missed some of your responses before. I welcome any practical references for repeated LMM effect size calculations. Note that I did not find Raudenbush or Snijders (as referenced in this technote http://www-01.ibm.com/support/docview.wss?uid=swg21474983) useful to that end. E From: Rudobeck, Emil (LLU)
Sent: Sunday, July 17, 2016 12:49 AM To: Ryan Black; [hidden email] Subject: RE: [BULK] [SPSSX-L] Mixed Model Effect Size Never mind my last message. I just stumbled upon this very same thread in the SPSSX online forum and saw responses to my message that I never received via email (I'm using only email
for this communication and hadn't checked out the forum much)...
Alright, so I reviewed the code more closely. I had mistakenly assumed it's a syntax for some complex calculation of an effect size (probably because I saw sigma and rho there), but all that syntax is just for generating a Gaussian dataset. Correct me if I'm wrong, but it just boils down to this: the residual variance of AR1 and ARMA11 structures is equivalent to the constant variance across the repeated measures (diagonal). Is that it? If so, it would be good to see some sort of a validating source on this because this method ignores other errors, such as the autocorrelation errors between the measurements (rho). Maybe that's not important, but I can't be sure without a reference. I'm also not clear why the syntax ends with two pseudo-r square calculations (equivalent output). You can only get one r square from this set of two models, but let me know if I'm missing something. I know that AR1 is not heterogeneous. In my original message, that was clipped, I was referring to ARH1. But we can get to that after clarifying the effect size for the homogeneous structures. The type of structure is dataset dependent and I have found a best fit with ARMA11, AR1, and even AD1 and TP structures. Emil From: Rudobeck, Emil (LLU)
Sent: Friday, July 15, 2016 8:30 PM To: Ryan Black; [hidden email] Subject: RE: [BULK] [SPSSX-L] Mixed Model Effect Size I think my last message was rejected because of attachments. I don't know what types of attachments are allowed, if any, but here is a link to the files:
https://drive.google.com/folderview?id=0B0AZLzY6tv15azBIaDU1bjZyU28&usp=sharing
And here is the message again: Hi Ryan, Thanks for the details. I am not that great with SPSS syntax. I think in your code you're just generating a dataset from random numbers and calculating based on those values, but I am not sure about the details, like why do you set x, e values = -99, unless it's just random designation. Also, rho = 0.5 presumably would be substituted with the actual rho value from unconditional the model. This example would be easier to understand if we get rid of the random number syntax. To that end, I have attached a simplified dataset along with the MIXED syntax for both AR1 and ARMA11. It would be helpful to see how your syntax for effect size would integrate with such a sample dataset to see exactly what is being calculated. Some questions: 1. A lot of times I encounter heterogeneous structure as the best fit. I'm assuming this effect size method applies not only to AR1, but to ARH1 as well. 2. My original issue was and still is with ARMA11. Are you aware of the methods for calculating the residual variance and/or the effect size for ARMA11 structures? In addition to rho, there is also phi for ARMA11, so AR1 effect size syntax probably can't be used for it. Thanks, Emil From: SPSSX(r) Discussion [[hidden email]] on behalf of Ryan Black [[hidden email]]
Sent: Wednesday, July 13, 2016 5:54 PM To: [hidden email] Subject: Re: [BULK] [SPSSX-L] Mixed Model Effect Size OP and others:
Although time is against me, I figured I'd estimate the pseudo R-squared estimate associated with the fixed effect of time from the model with an AR1 specification I posted earlier today.
For this model, pseudo R-squared can be calculated as:
(unconditional residual variance - conditional residual variance) / unconditional residual variance
OR
1 - (conditional residual variance / unconditional residual variance)
Below is the code.
Feel free to write back if you have questions. I'll try to respond as time permits.
HTH.
Ryan
--
*Generate Data for Mixed Model with AR1 specification.
set seed 65923454. new file.
inp pro. compute subject=-99.
compute time = -99. compute x1 = -99. compute x2 = -99. compute x3 = -99. compute e1 = -99. compute e2 = -99. compute e3 = -99. compute sigma = 1. compute rho = 0.50. compute a11 = 1. compute a21 = rho. compute a31 = rho**2. compute a22 = sqrt(1 - rho**2). compute a32 = rho*sqrt(1 - rho**2). compute a33 = sqrt(1 - rho**2). leave subject to a33.
loop subject= 1 to 1000.
compute x1 = rv.normal(0,1).
compute x2 = rv.normal(0,1). compute x3 = rv.normal(0,1). compute e1 = sigma * a11*x1.
compute e2 = sigma * (a21*x1 + a22*x2). compute e3 = sigma * (a31*x1 + a32*x2 + a33*x3). loop time = 1 to 3.
compute y = 1.5 + 1.2*(time=1) + 0.9*(time=2) + e1*(time=1) +
e2*(time=2) + e3*(time=3). end case. end loop. end loop. end file. end inp pro. exe. delete variables x1 x2 x3 sigma rho a11 a21 a31 a22 a32 a33 e1 e2 e3.
MIXED y BY time
/FIXED= | SSTYPE(3) /METHOD=REML /PRINT=R SOLUTION /REPEATED=time | SUBJECT(subject) COVTYPE(AR1). MIXED y BY time
/FIXED=time | SSTYPE(3) /METHOD=REML /PRINT=R SOLUTION /REPEATED=time | SUBJECT(subject) COVTYPE(AR1). compute pseudo_R_squared_1 = (1.319641 - 1.035615) / 1.319641.
compute pseudo_R_squared_2 = 1 - (1.035615 / 1.319641). execute. On Wed, Jul 13, 2016 at 1:49 PM, Ryan Black
<ryan.andrew.black@...> wrote:
WARNING: Please be vigilant when opening emails that appear to be the least bit out of the ordinary, e.g. someone you usually don’t hear from, or attachments you usually don’t receive or didn’t expect, requests to click links or log into systems, etc. If you receive suspicious emails, please do not open attachments or links and immediately forward the suspicious email to [hidden email] and then delete the suspicious email. |
Maybe Jack did; see:
http://journal.frontiersin.org/article/10.3389/fpsyg.2012.00111/full -Mike Palij New York University [hidden email] P.S. Regarding AuC/A' measures, context is important as are assumptions. In the traditional Signal Detection Theory model, nobody used AuC/A' because if the assumption of normal distributions with homogeneous variances was met, then the preferred measure of "sensitivity" or correct performance was d-prime or d' which many people will recognize as being similar to Cohen's d, that is, d = (mu1 - mu2)/sigma. But old school SDT typically would set up different experiemental conditions where the criterion Beta (a threshold value on a likelihood ratio continuum) would vary, produce a different hit rate and false alarm rate for each condition. Each pair of values would fall on the isosentitivity curve (constant sensitive or d') or the ROC curve. After several thousand trials in a psychophysics experiment, assumptions were more or less met, and SDT gave a good description of what was going on. However, when it was learned that variances would be heterogeneous and that bias might be related to sensitivity, other measures were developed. Pertinent to the current discussion is the situation where one only has one pair of hit and false alarm rates -- this doesn't really allow one to have a classic isosensotivity curve and d-prime is not appropriate. This situation would occur often in SDT analysis of recognition memory and "nonparametric" sensitivity were developed to replace d-prime, among them A-prime or A' or AuC (similar new statistics for the crieterion beta were developed such as B'). Research in psychophysics and recognition memory could use either a binary response (e.g., Yes-No) or rating scale response (e.g., a 5 point rating scale) where each value had an associated hit rate and false alarm rate -- this allowed one to determine the isosensitivity curve from one person in one experimental condition. SDT has continued to be developed (Larry DeCarlo at Teacher's College at Columbia has been some work in this area, such as linking SDT to the general linear model and the use of latent variables). John Swets brought SDT analysis to radiology in order to provide quantitative methods for assessing radiologist's ability to correctly read X-rays and other scans, and other applications spread as the researcher became familiar with SDT concepts and analysis. Diagnosis in clinical psycholog using SDT has been suggested by variious people such as: McFall, R. M., & Treat, T. A. (1999). Quantifying the information value of clinical assessments with signal detection theory. Annual review of psychology, 50(1), 215-241. But clinicians tend to be people people instead of math people though Jack Cohen got his Ph.D. in clinical psych and Peter Bentle of SEM and EQS fame was another clinical psychologist gone math. So, the real question comes down to whether a SDT type analysis with AuC or other measures of sensitivity or correct performance meet certain assumptions or whether one modifies the nature of the analysis with different assumptions (in grad school I took a course with a math psychologist who presented SDT not with normal distributions but with exponential distributions -- even though I was familiar with traditional SDT I, like everyone in the course, was completely lost). Both ROC and Kappa and related measures have become popular in AI as measurs of performance in machine learning, pattern recognition, and so on, hence the interest in multiclass analysis. Quotomg from the entry on ROC analysis in the Handbook of Machine Learning: |Future Directions | |ROC analysis in its original form is restricted to binary classification, |and its extension to more than two classes gives rise to many open |problems. c-class ROC analysis requires c(c1) dimensions, in |order to distinguish each possible misclassification type. Srinivasan |proved that basic concepts such as the ROC polytope and its linearly |interpolated convex hull generalise to the c-class case [15]. |In theory, the volume under the ROC polytope can be employed for |assessing the quality of a multi-class classifier [7], but this volume is |hard to compute as – unlike the two-class case, where the segments |of an ROC curve can simply be enumerated in O(nlogn) time by |sorting the n examples on their score [5, 9] – there is no simple way |to enumerate the ROC polytope. Mossman considers the special |case of 3-class ROC analysis, where for each class the two possible |misclassifications are treated equally (a so-called one-versus-rest |scenario) [13]. Hand and Till propose the average of all one-versus-rest |AUCs as an approximation of the area under the ROC polytope [11]. |Various algorithms for minimising a classifier’s misclassification costs |by reweighting the classes are considered in [12, 1]. From: Flach, P. A. (2011). ROC analysis. In Encyclopedia of machine learning (pp. 869-875). Springer US. The case of 3 Even the weather folks are getting involved; see: Wandishin, M. S., & Mullen, S. J. (2009). Multiclass ROC analysis. Weather and Forecasting, 24(2), 530-547. http://journals.ametsoc.org/doi/abs/10.1175/2008WAF2222119.1 ===================== To manage your subscription to SPSSX-L, send a message to [hidden email] (not to SPSSX-L), with no body text except the command. To leave the list, send the command SIGNOFF SPSSX-L For a list of commands to manage subscriptions, send the command INFO REFCARD |
In reply to this post by Rudobeck, Emil (LLU)
Emil: You've said quite a bit in your most recent post, some of which I disagree. Let me focus in on your question about estimating Cohen's D. Suppose you wanted calculate Cohen's d as a standardized estimate of the difference in means between two groups in a linear mixed model (fixed effect component includes the group effect only) that assumes a first-order autoregressive residual correlation matrix. One reasonable solution would be: d = mean diff / sd where, mean diff = estimated mean difference between groups i and j sd = square root of any one of the elements on the diagonal of the residual correlation matrix In more complex cases, where there are unequal residual diagonal values for example, the story does become a bit more complicated. Broader point. Calculating an effect size in a linear mixed model depends on many factors such as the type of effect size you are interested in and for which predictor(s), as well as the random effects and residual variance covariance matrices. For example, if one employed a linear mixed model with two levels, it would be possible to calculate something analogous to eta-squared related to a level 2 predictor or level 1 predictor. If the responses you are receiving are not to your satisfaction, then perhaps you might consider doing some of the heavy lifting to arrive at an answer that is satisfactory to you, and then submit your work for review by the scientific community. Ryan On Mon, Nov 28, 2016 at 11:55 PM, Rudobeck, Emil (LLU) <[hidden email]> wrote:
|
Ryan,
=====================
To manage your subscription to SPSSX-L, send a message to
[hidden email] (not to SPSSX-L), with no body text except the
command. To leave the list, send the command
SIGNOFF SPSSX-L
For a list of commands to manage subscriptions, send the command
INFO REFCARD
Perhaps my broad approach to finding out the best measure of effect size is inefficient given that it's an area of active research. Calculating Cohen's d makes intuitive sense since there are other methods similar to what you mentioned and I think for an H structure I can just stick to a timepoint. Maybe we could focus on the overall variance of a factor based on published work. Equation 26 in the Nakagawa paper (link below) refers to r squared as a ratio of several variances: 1. Var(f) = Seems to be the variance obtained from the R matrix when the model is run without any random effects 2. Var(y) = Group-specific variance 3. Var(a) = Between-group variance 4. Var(e) = Within-group variance Ignoring the issue of heterogeneous covariances, do you know how these variances can be obtained from the SPSS MIXED output? When Repeated isn't specified, I would think the variances would be similar to methods when calculating ICC (within error = residual variance, between error = intercept variance). But once a structure like ARMA11 is specified, it's unclear. Maybe because the method isn't applicable anymore? Reference: Nakagawa 2013: http://onlinelibrary.wiley.com/doi/10.1111/j.2041-210x.2012.00261.x/full From: SPSSX(r) Discussion [[hidden email]] on behalf of Ryan Black [[hidden email]]
Sent: Tuesday, November 29, 2016 5:35 AM To: [hidden email] Subject: Re: [BULK] [SPSSX-L] Mixed Model Effect Size Emil:
You've said quite a bit in your most recent post, some of which I disagree.
Let me focus in on your question about estimating Cohen's D. Suppose you wanted calculate Cohen's d as a standardized estimate of the difference in means between two groups in a linear mixed model (fixed effect component includes the group effect only)
that assumes a first-order autoregressive residual correlation matrix.
One reasonable solution would be:
d = mean diff / sd
where,
mean diff = estimated mean difference between groups i and j
sd = square root of any one of the elements on the diagonal of the residual correlation matrix
In more complex cases, where there are unequal residual diagonal values for example, the story does become a bit more complicated.
Broader point. Calculating an effect size in a linear mixed model depends on many factors such as the type of effect size you are interested in and for which predictor(s), as well as the random effects and residual variance covariance matrices. For example,
if one employed a linear mixed model with two levels, it would be possible to calculate something analogous to eta-squared related to a level 2 predictor or level 1 predictor.
If the responses you are receiving are not to your satisfaction, then perhaps you might consider doing some of the heavy lifting to arrive at an answer that is satisfactory to you, and then submit your work for review by the scientific community.
Ryan
On Mon, Nov 28, 2016 at 11:55 PM, Rudobeck, Emil (LLU)
<erudobeck@...> wrote:
|
Emil: That formula assumes a nested random effects linear mixed model. It's certainly possible to obtain those estimated variance components of the random factors using the MIXED procedure in SPSS if you construct the model with nested random effects. They also show explicitly how to obtain the fixed effects variance. I can show how at a later date but I'm swamped for a while. Ryan Sent from my iPhone ===================== To manage your subscription to SPSSX-L, send a message to [hidden email] (not to SPSSX-L), with no body text except the command. To leave the list, send the command SIGNOFF SPSSX-L For a list of commands to manage subscriptions, send the command INFO REFCARD |
Free forum by Nabble | Edit this page |