Hello, I have two treatments (treated and control) with 20 animals each treatment.
And I have a reactivity score (0 no reactive and 1reactive) from 10 days (day 1, 3, 5, 7, ..., 17, 19).
I would like to compare A and B by day just for escore 1.
So the objective is see if the treated group had less reactivity than group control on days 1 till 19.
They told me to use proc genmod or glimmix.
I have measure repeated in time: I have 10 information (1 information by day collected) per animal.
Can somebody help me, please?
proc genmod data=fake descending
order=data;
class treatment animal;
model reactivity=treatment day/ d=bin link=logit;
repeated subject=animal / type=cs;
estimate "treatment" treatment 1 1 / exp;
estimate "day" day / exp;
run;
The above code fits a GEE model with the binomial distribution.
Thank you,
Rajesh.
proc genmod data=fake descending
order=data;
class treatment animal;
model reactivity=treatment day/ d=bin link=logit;
repeated subject=animal / type=cs;
estimate "treatment" treatment 1 1 / exp;
estimate "day" day / exp;
run;
The above code fits a GEE model with the binomial distribution.
Thank you,
Rajesh.
Hi Rajesh, Thank you for your fast reply! Very kind of you, thank you very much for helping me.
So I tried to do it and it works for few, but I got errors in others.
I put it like this:
proc genmod data=C descending order=data;
class group cow dl po;
model react=group dl po/ d=bin link=logit;
repeated subject=cow/ type=cs;
estimate "group" group1 1 / exp;
estimate "dl" dl1 1/ exp;
run;
group is control and treat
dl is the day of lactation (1, 3, 5, 7, 9, 11, 13, 15, 17 and 19)
po is parturition order
react is the reactivity scored as 0 (without reactivity) and 1 (with reactivity).
In my table, I will just compare the score 1 (with reactivity). As example below:
DL  Reactivity  
Group Treated  Group Control  P value  
01  50  80 

03  30  65 

05  20  35 

I would like to have the P value.
Sometimes the difference between groups is large and the P value isn't significative.
And sometimes have errors like this:
NOTE: PROC GENMOD is modeling the probability that react='1'.
WARNING: The negative of the Hessian is not positive definite. The convergence is questionable.
WARNING: The procedure is continuing but the validity of the model fit is questionable.
WARNING: The specified model did not converge.
WARNING: Negative of Hessian not positive definite.
WARNING: The generalized Hessian matrix is not positive definite. Iteration will be terminated.
ERROR: Error in parameter estimate covariance computation.
ERROR: Error in estimation routine.
I'm bit new in this... and a desesperate student.
Please, if you can help me, I am really glad.
Thank you again.
Aska.
proc genmod data=C descending order=data;
class group cow po;
model react=group dl po/ d=bin link=logit;
repeated subject=cow/ type=cs;
estimate "group" group1 1 / exp;
estimate "dl" dl/ exp;
run;
The 'dl' variable need not be mentioned in the class statement. If you do not mention it in the class statement, it will be treated as a continuous variable. Also, if you do not need to control for 'po' variable, please remove it from the model. Please check if the model converges with the above code.
Based on the error that you are getting, you have specified a more complex model than your data can support.
Hope this helps!
Thank you,
Rajesh.
Thank you for you reply and help!
It works now, thank you for you explanation.
But, abusing again your goodwill and knowledge, now the results havent the information of estimate, SE, Z and Pr>Z values for each DL...
Thank you again!
Aska.
It will not give estimates for each time point separately. Now the interpretation becomes "change in reactivity for a day increase in time". If you want to interpret it for 2 days increase multiply the estimate for 'dl' by 2 and exponentiate it to get the appropriate change. Else, replace the relevant estimate statement with the below code.
estimate 'dl' dl 2/exp;
Thanks,
Rajesh.
Hello Rajesh, thank you.
This is the last thing, I promise. When I put: estimate "dl" dl 1 1/ exp;
What this "1 1" means?
Thank you!
Aska.
If you included 'dl' in the class statement, the below SAS code will apply.
estimate "dl" dl 1 1/ exp;
It means for a change in day from 1 to 3/1 to 5 (this will change based on the reference chosen but the concept is this), the change (increase/decrease) in the odds of having reactivity=1.
If you did not include 'dl' in the class statement, the estimate statement should be:
estimate "dl" dl 1/ exp;
It means, on average, for a unit increase in day, the change (increase/decrease) in the odds of having reactivity=1.
Thank you,
Rajesh.
Thank you Rajesh again.
So that mean that is comparing day with day by genmod procedure.
And probably I have to change my PROC because what I would like to do is compare the difference between groups for each day and not compare between days for presence of reactivity.
Thank you for helping me!!! It really helped me!!
Aska
Build your skills. Make connections. Enjoy creative freedom. Maybe change the world. Registration is now open through August 30th. Visit the SAS Hackathon homepage.
Register today!SAS' Charu Shankar shares her PROC SQL expertise by showing you how to master the WHERE clause using real winter weather data.
Find more tutorials on the SAS Users YouTube channel.