Back to the Top
hey all,
This is a question about how to calculate the standard error of the adjusted treatment means when
it comes to two covariates and two way analysis of covariance case. It might not related with pkpd, but
I hope someone could help me, since there might be someone who is good at mathematics. I spent a lot of time on it.
This question focuses on linear algebra and arithmetic. That means I do not want to use any software.
the SAS code list below and following 4 questions :
this is 2*4 case
treatment=2
center=4
--
two covariates: AGE BPDIA0
= SAS CODE
--
/*Example 11.2 (ANCOVA in MultiCenter Hypertension Study) */
DATA BP;
INPUT TRT $ PAT AGE BPDIA0 BPDIACH .-at-..aaa.;
CENTER = INT(PAT/100);
DATALINES;
/*Example 11.2 (ANCOVA in MultiCenter Hypertension Study) */
DATA BP;
INPUT TRT $ PAT AGE BPDIA0 BPDIACH .-a-..at.;
CENTER = INT(PAT/100);
DATALINES;
SOLE 101 55 102 -6 SOLE 103 68 115 -10
SOLE 104 45 97 -2 SOLE 106 69 107 -5
SOLE 109 54 115 -7 SOLE 202 58 99 -4
SOLE 203 62 119 -11 SOLE 204 51 107 -9
SOLE 206 47 96 0 SOLE 208 61 98 -6
SOLE 211 40 110 5 SOLE 212 36 103 -3
SOLE 214 58 109 10 SOLE 216 64 119 -12
SOLE 217 55 104 4 SOLE 219 54 97 -5
SOLE 223 39 95 -8 SOLE 301 49 115 -3
SOLE 302 46 105 4 SOLE 304 59 116 -10
SOLE 307 42 108 -5 SOLE 308 65 101 -7
SOLE 311 68 102 -5 SOLE 313 57 110 -8
SOLE 402 73 119 -12 SOLE 404 48 99 -7
SOLE 406 53 117 0 SOLE 407 46 96 4
SOLE 409 60 118 -15 SOLE 412 66 104 -3
SOLE 414 59 115 -4 SOLE 415 41 109 3
SOLE 418 53 116 -10 SOLE 421 57 100 -8
SOLE 424 52 103 0 SOLE 425 41 95 9
COMB 102 68 112 -10 COMB 105 64 105 -4
COMB 107 48 107 -9 COMB 108 60 107 -5
COMB 201 44 109 -7 COMB 205 53 99 0
COMB 207 48 107 -9 COMB 209 65 106 -3
COMB 210 79 108 -6 COMB 213 45 110 -9
COMB 215 44 93 0 COMB 218 62 99 -5
COMB 220 59 119 -14 COMB 221 50 104 -6
COMB 222 61 107 -14 COMB 224 36 95 3
COMB 225 34 95 -9 COMB 303 49 115 -8
COMB 305 57 115 7 COMB 306 58 116 -20
COMB 309 43 108 -7 COMB 310 66 101 -5
COMB 312 55 102 -9 COMB 314 70 110 -16
COMB 315 52 104 -8 COMB 401 42 119 -15
COMB 403 49 109 -6 COMB 405 53 117 4
COMB 408 55 96 6 COMB 410 65 98 -8
COMB 411 75 108 -7 COMB 413 59 104 -10
COMB 416 60 115 -16 COMB 417 53 109 3
COMB 419 43 96 -14 COMB 420 77 100 -5
COMB 422 55 103 -8 COMB 423 66 115 -14
COMB 426 38 105 -3
;
PROC PRINT DATA = BP;
TITLE1 'Analysis of Covariance';
TITLE2 'Example 11.2: ANCOVA in Multi-Center Hypertension Study';
RUN;
PROC SORT DATA = BP;
BY TRT CENTER;
PROC MEANS MEAN STD N MIN MAX DATA = BP;
BY TRT;
VAR AGE BPDIA0 BPDIACH;
RUN;
PROC GLM DATA = BP;
CLASSES TRT CENTER;
MODEL BPDIACH = TRT CENTER AGE BPDIA0 / SOLUTION;
LSMEANS TRT CENTER/ STDERR PDIFF;
TITLE3 'GLM Model Including Effects: TRT, CENTER, AGE, BPDIA0';
QUIT;
RUN;
= --
SAS OUTPUT
--
The MEANS Procedure
Variable Mean Std Dev N Minimum Maximum
--
AGE 55.3846154 11.0706235 39 34.0000000 79.0000000
BPDIA0 106.3333333 7.0012530 39 93.0000000 119.0000000
BPDIACH -6.8205128 6.2318021 39 -20.0000000 7.0000000
--
-- TRT=SOLE --
Variable Mean Std Dev N Minimum Maximum
--
AGE 54.1944444 9.4923110 36 36.0000000 73.0000000
BPDIA0 106.6666667 7.9677923 36 95.0000000 119.0000000
BPDIACH -4.0555556 5.9997354 36 -15.0000000 10.0000000
--
--
The GLM Procedure
Class Level Information
Class Levels Values
TRT 2 COMB SOLE
CENTER 4 1 2 3 4
Number of Observations Read 75
Number of Observations Used 75
The GLM Procedure
Least Squares Means
H0:LSMean1 BPDIACH Standard H0:LSMEAN=0 LSMean2
TRT LSMEAN Error Pr > |t| Pr > |t|
COMB -6.93129242 0.98013179 <.0001 0.0435
SOLE -4.19121438 0.99307483 <.0001
BPDIACH Standard LSMEAN
CENTER LSMEAN Error Pr > |t| Number
1 -5.86408602 1.93424320 0.0034 1
2 -5.55003909 1.17247036 <.0001 2
3 -5.93099721 1.49367504 0.0002 3
4 -4.89989128 1.12770583 <.0001 4
--
= I'd found there are book named, the author's offered
the variance of adjusted treatment means at page 300
the formula is:
Var(Y_adj)=(T'X_c)(X'X)^- (T'X_c)' Sigma^2
T'=n*(a+b)=(2+4)*75 matrix
X_c=75*[(2+4)+2] matrix
(T'X_c) should be a 6*8 matrix
X=(a+b)+two covariate =[(2+4)+2]*n
so (X'X) is a 8*8 matrix
(T'X_c)(X'X)^-(T'X_c)'
should be 6*6 matrix
but i did not get the same result after computing my VBA code.
perhaps there might be some misunderstanding about the above formula.
--
Question1: what is b?
--page 297
T=1/b*|matrix|
matrix is list below
1 0 ... 0
.. . ... .
1 . ... 0
0 1 ... 0
.. . .
0 1 ... 0
.... ... ...
0 0 ... 1
.. . .. .
0 0 ... 1
b is observations per treatment.
in above case, COMB has 39 observations, SOLE has 36 observations.
so b=?
and there also have 4 levels of CENTER
what b should be?
--
Question2: T matrix
--
the matrix in question 1
is a*n matrix, n is observations, a is levels of class variable. but for two way and two covariates analysis
is T matrix like this?
(a+b)*n =(2+4)*75 matrix: |.....a....| |....b.......|
1 0 ... 0 1 0 ... 0
.. . ... . . . ... .
1 . ... 0 1 . ... 0
0 1 ... 0 0 1 ... 0
.. . .
0 1 ... 0 0 1 ... 0
.... ... ... ... ... ...
0 0 ... 1 0 0 ... 1
.. . .. . . . .. .
0 0 ... 1 0 0 ... 1
a=2 is levels of TRT, b=4 is levels of CENTER.
--
Question3: X_c matrix
--
in page 300, equation 9.46
Var(Y_adj)=(T'X_c)(X'X)^- (T'X_c)' Sigma^2
if no class variable X 4*n matrix
group center age bpdia0 1 1 55 102 1 1 45 97 .... ... ... ... ...
2 1 68 112 2 1 48 107 .... ... ... ... ... 1 3 46 105 .... ... ... ... ...
1 4 73 119 2 1 68 112 2 2 44 109
.... ... ... ... 2 3 52 104 .... ... ... ...
2 4 49 109 2 4 38 105 .... ... ... ...
if HAVE two class variables X[(2+4)+2]*75 matrix
--
|a||...b....|
1 0 1 0 0 0 55 102 1 0 1 0 0 0 45 97 .... ... ... ... ...
0 1 1 0 0 0 68 112 0 1 1 0 0 0 48 107 .... ... ... ... ... 1 0 0 0 1 0 46 105 .... ... ... ... ...
1 0 0 1 0 0 73 119 0 1 1 0 0 0 68 112 0 1 0 1 0 0 44 109
.... ... ... ... 0 1 0 0 1 0 52 104 .... ... ... ...
0 1 0 0 0 1 49 109 0 1 0 0 0 1 38 105 .... ... ... ...
--
is X_c [(2+4)+2] *75 matrix (two covariates=0 )
|a||....b...|
1 0 1 0 0 0 0 0 1 0 1 0 0 0 0 0 .... ... ... ... ...
0 1 1 0 0 0 0 0 0 1 1 0 0 0 0 0
.... ... ... ... ... 1 0 0 0 1 0 0 0 .... ... ... ... ...
1 0 0 1 0 0 0 0
0 1 1 0 0 0 0 0
0 1 0 1 0 0 0 0
.... ... ... ... 0 1 0 0 1 0 0 0
.... ... ... ...
0 1 0 0 0 1 0 0 0 1 0 0 0 1 0 0
.... ... ... ...
is it right?
--
Question 4: the classical formula
--
In bookpage 300
Exx=sigma i f= 1 to a, sigma j=1 to b [X_ij - mean(X_i)-mean(X_.j)+mean(X_..)]
for my understanding
mean(X..)=1/2(mean(X_i.)+mean(X_.j))
so, Exx=sigma i = 1 to a, sigma j=1 to b [X_ij -mean(X_..)]
--
isn't it?
sorry about so long an email.
Regards
Kun Wang Ph.D
Center for Drug Clinical Research
Shanghai University of Chinese Medicine
1200 Cailun Road, 201203 Shanghai, China
Tel: 86-21-51322420
Mb: 13816298627
Back to the Top
The following message was posted to: PharmPK
I can't answer the question directly, but if you don't get a good answer
from this illustrious gropup, you might like to try posting to
StackOverflow as a programming related question
http://stackoverflow.com/questions/tagged/sas
or the new Statistics
StackExchange beta Q&A site
http://stats.stackexchange.com/questions
as a pure stats question.
Regards,
Paul.
--
Paul Hurley
Want to post a follow-up message on this topic?
If this link does not work with your browser send a follow-up message to PharmPK@boomer.org with "Standard error of adjusted treatment means" as the subject | Support PharmPK by using the |
Copyright 1995-2011 David W. A. Bourne (david@boomer.org)