Date: Wed, 18 Apr 2007 16:20:01 -0700
Reply-To: Dale McLerran <stringplayer_2@YAHOO.COM>
Sender: "SAS(r) Discussion" <SAS-L@LISTSERV.UGA.EDU>
From: Dale McLerran <stringplayer_2@YAHOO.COM>
Subject: Re: PROC IML
In-Reply-To: <113678.50540.qm@web63503.mail.re1.yahoo.com>
Content-Type: text/plain; charset=iso-8859-1
--- Yuanhui Zhu <yzhu97@yahoo.com> wrote:
> Hi Dale,
> I have a multivariate response with irregularly spaced
> time points which is unique for every individual. So, I want to
> explicitly perform the MLE. Any thoughts ?
>
> Thanks,
> Y. Zhu
>
>
How about the following:
proc mixed data=mydata;
class response subject;
model y = response response*x1 response*x2 ... / noint s;
repeated / subject=subject type=sp(pow)(time) group=response;
random response / subject=subject type=un;
run;
where RESPONSE is an indicator of which response variable is
being modeled. So, if you originally have a multivariate response
Y1, Y2, Y3
Subject Y1 Y2 Y3 time
1 2.87 9.82 7.05 1
1 3.19 8.85 7.34 16
1 ... ... ... ...
2 4.23 12.06 9.81 1
2 6.12 11.40 9.21 23
2 ... ... ... ...
... ... ... ... ...
then you would reshape the data to appear as
Subject Y Response Time
1 2.87 1 1
1 9.82 2 1
1 7.05 3 1
1 3.19 1 16
1 8.85 2 16
1 7.34 3 16
2 4.23 1 1
2 12.06 2 1
2 9.81 3 1
2 6.12 1 23
2 11.40 2 23
2 9.21 3 23
...
So then, what is done by the above code? Let's start with the MODEL
statement where we have removed the intercept and then specified
the interaction of RESPONSE with each predictor variable. With
the intercept removed and RESPONSE specified on the CLASS statement,
then the dummy variable expansion of RESPONSE will provide us with
intercept estimates for each response variable. The interaction
terms will provide separate slope estimates on each response for
each of the predictor variables X1, X2, ....
Now, what about the REPEATED statement. We use a spatial
covariance structure which allows separate covariance structures
for each of the response variables. Note that if we fit the model
proc mixed data=mydata;
class response subject;
model y = response response*x1 response*x2 ... / noint s;
repeated / subject=subject type=sp(pow)(time) group=response;
run;
the parameter estimates should be identical to estimates obtained
for
proc mixed data=mydata;
by response;
class subject;
model y = x1 x2 ... / s;
repeated / subject=subject type=sp(pow)(time);
run;
In this last model, we make no attempt to treat the response as
a multivariate vector. We fit separated models for each response
variable, and for each response we allow a spatial covariance
structure in order to account for the different data collection
times for each individual.
So, that brings us to the final statement in my grand model, the
RANDOM statement. The random statement enables us to account for
covariance in the multivariate response measurements within
subjects. Observe in the (limited) data for our first two subjects
that, conditional on the response variable, the values of the
response in subject 1 are all less than the values of the response
in subject 2. Each subject has their own expectation. The
RANDOM statement allows us to account for subject-specific
expectations on each of the response variables. The subject-
specific expectations would modify the covariance between the
responses.
If each subject had responses at only two time points and our model
has just two response variables, then the above model would give us
a covariance structure as developed below:
Let Yit be response i at time t.
Y11 ~ N( n11 + g1, | s1^2 (s1^2)*rho1^(t2-t1) | )
Y12 ( n12 + g1, | (s1^2)*rho1^(t2-t1) s1^2 | )
Y21 ~ N( n21 + g2, | s2^2 (s2^2)*rho2^(t2-t1) | )
Y22 ( n22 + g2, | (s2^2)*rho2^(t2-t1) s2^2 | )
g1 ~ N( 0, | S11 S12 | )
g2 ( 0, | S12 S22 | )
nit = etaij is the fixed effect for response i at time t
g1 is the subject-specific random effect for response 1
g2 is the subject-specific random effect for response 2
We can also write Yit = nit + gi + eit were nit is the fixed
effect for response i at time t, gi is the random effect for
response i, and eit is the residual error for response i at
time t. From the above, we find that the covariance structure
for the model I present is
V(Y11) = V(Y12) = s1^2 + S11
V(Y21) = V(Y22) = s2^2 + S22
Diagonal variance is variance of residual error plus variance
of random effect. For each variable, the residual error
variance is identical over time and we have the same random
effect so within variable over time variance structure is
constant on the diagonal.
cov(Y11, Y12) = (s1^2) * rho1^(t2-t1) + S11
cov(Y21, Y22) = (s2^2) * rho2^(t2-t1) + S22
Within response covariances over time - cov through repeated
measures and shared random effect
cov(Y11, Y21) = S12
cov(Y11, Y22) = S12
cov(Y12, Y21) = S12
cov(Y12, Y22) = S12
Between response covariances - covariance through covariance
of g1 and g2 random effects only
If you want to estimate something like a UN@sp(pow)(time) repeated
measures model where you assume that the spatial power function
is the same for each response but the variance is different for
each of the responses so that you would have
V(Y11) = V(Y12) = S11
V(Y21) = V(Y22) = S22
cov(Y11, Y12) = S11 * rho^(t2-t1)
cov(Y21, Y22) = S22 * rho^(t2-t1)
cov(Y11, Y21) = S12
cov(Y11, Y22) = S12 * rho(t2-t1)
cov(Y12, Y21) = S12 * rho(t2-t1)
cov(Y12, Y22) = S12
then you will need to write IML code.
Dale
---------------------------------------
Dale McLerran
Fred Hutchinson Cancer Research Center
mailto: dmclerra@NO_SPAMfhcrc.org
Ph: (206) 667-2926
Fax: (206) 667-5977
---------------------------------------
__________________________________________________
Do You Yahoo!?
Tired of spam? Yahoo! Mail has the best spam protection around
http://mail.yahoo.com
