ABBREVIATIONS:

A Active Drug

B Placebo

Baseline BP B0

BP Blood Pressure [in mmHG]

CV Cardio Vascular Disease

DOD Difference of Difference (T test)

LM linear model

MDS Multi-dimension sclaing

TRT Treatment either placebo or Drug Beta blocker

RCT Randomized Clinical Trial

RM Repeated Measures

RTM Regression to the Mean

SD Standard deviation

| Given that

:: Packages R

GOAL

The study of longitudinal data is a challenge for Scientists ,Students and practitioner because one wants to make inference and modelling but often failed to account for RM in their classic statistical tests : Assumptions are broken and is a flaw in publications : As result coefficient might be still reflect true value but the SE of coefficient or test or inflated¹.

The main subject here is how including a baseline as a co variate act on your “Linear parameters” model prediction.

We present a short analysis as follow:

RCT & Data. Wide vs Long format
Graphing longitudinal data [ID trajectories]²
Some word on variability in time
Guideline basic-naive linear modelling: T-Test on endpoint 4 months and LM
Linear model (with/without baseline inclusion)
Linear mixed model

RCT and DATA

Beta blocker and Reserpine like drugs³ has been proven significantly efficient in CV diseases (Stroke,Infarction,BP reduction etc), this since the 60’s which was the golden age of these discoveries.

The physiological mechanisms are:

B1 or VMAST calcic antagonist blocking agent⁴ reduce contractile ventricular pumping hence reducing the load on arterial and peripheral pipe systems
Less blood flow
Inhibition of renine-angiotensine II at kidney receptor site -An reduce autonmous nervous system activity.
Vaso-dilatation of vessels and much more effects

Of course multi-targeted spots in our cells receptors (via Calcium or B blockade chained) result in ends of a major drop in our BP.

The data presented here come from a study: Veterans Administration Cooperative Study Group on Anti hypertensive Agents (VA Study Group (1967, D. Freis):
Endpoints Y_it: Diastolic BP Supine horizontal in [mmHg]

Beside our BP follwows cycles and might be influenced by: Natural Rhythm (“random”) , Emotions-Stress (“random”) , Cafeine-Alcholl-drugs-exercises consumption (can be controlled) among others.
Randomization between placebo and TRT were done.
Questions or interest : Does the new Drug reduce the Diastolic BP with time?

Covariates : Sex,Age and a Baseline [BP0] record.

Response: Y_it are the BP at 1 months intervals: [BP0] BP1 BP2 BP3 BP4

Duration 4 months follow-up from Baseline inwards.

Note: Data from ID1 was removed to create unbalanced design and show up difficulties in running within subject ANOVA..

Note that Baseline is one a covariate once as a real measure in time=0 and we will explore how this possibilities act for inferences and modelling.

Know your Data: Systematic review

The data come form the book CLINICAL DATA ANALYSIS USING R (Chen & Peace , 2011) data extracted from the 1967 Trial: Some question arouse:
Does the author have copied partially the data? Does it exist an metadata of the DOE ?

Such questions must be answered before proceeding to the analysis.

Note that 20 id in Placebo group B and 20[-1]⁵ id in trt group [A] seems to be a phase II trial as n is low and the reader should bear in mind at this stage of analysis non significant result might come from the too low power of the study design, especially in longitudinal data.

In present study reader must be aware that regression to the mean might occur when :

1: When part of a population with abnormal value are selected

2: Random process occur under Normal law

3: In repeated measures RTM is exclusively noticeable

Note : could be attributed to the Placebo effect from practitioner

By RCT design we expect that RTM is the same magnitude in both arms and hence is accounted for in assignment of trt.

From the author:

We won’t account for sake of simplicity in the design, study of within between subject and LMM will be exposed in other blog Using lm in Longitudinal data is a flaw as RM are correlated but coefficients are usually usable but the SE for inferences.

We also try to use different plotting packages i.e ggplot2 and lattice to present the reader the pro and cons of each one.(let at your appraisal)

Data summary

## [1] "GROUP A= TRTEATMENT B= PLACEBO"

RAW DATA BP
id	trt	baseline	bp1	bp2	bp3	bp4	age	sex
2	Drug	116	113	112	103	101	51	M
3	Drug	119	115	113	104	98	48	F
4	Drug	115	113	112	109	101	42	F
5	Drug	116	112	107	104	105	49	M
6	Drug	117	112	113	104	102	47	M
7	Drug	118	111	100	109	99	50	F
8	Drug	120	115	113	102	102	61	M
9	Drug	114	112	113	109	103	43	M
10	Drug	115	113	108	106	97	51	M
11	Drug	117	112	110	109	101	47	F
12	Drug	116	115	113	109	102	45	M
13	Drug	119	117	110	106	104	54	F
14	Drug	118	115	113	102	99	52	M
15	Drug	115	112	108	105	102	42	M
16	Drug	114	111	111	107	100	44	F
17	Drug	117	114	110	108	102	48	M
18	Drug	120	115	113	107	103	63	F
19	Drug	114	113	109	104	100	41	M
20	Drug	117	115	113	109	101	51	M
21	Placebo	114	115	113	111	113	39	M
22	Placebo	116	114	114	109	110	40	F
23	Placebo	114	115	113	111	109	39	F
24	Placebo	114	115	113	114	115	38	M
25	Placebo	116	113	113	109	109	39	F
26	Placebo	114	115	114	111	110	41	M
27	Placebo	119	118	118	117	115	56	F
28	Placebo	118	117	117	116	112	56	M
29	Placebo	114	113	113	109	108	38	M
30	Placebo	120	115	113	113	113	57	M
31	Placebo	117	115	113	114	115	47	F
32	Placebo	118	114	112	109	110	48	M
33	Placebo	121	119	117	114	115	61	F
34	Placebo	116	115	116	114	111	49	M
35	Placebo	118	118	113	113	112	52	M
36	Placebo	119	115	115	114	111	55	F
37	Placebo	116	114	113	109	109	45	F
38	Placebo	116	115	114	114	112	42	M
39	Placebo	117	115	113	114	115	49	F
40	Placebo	118	114	114	114	115	50	F

SUMMARY WIDE DATA BP.CSV
id	trt	baseline	bp1	bp2	bp3	bp4	age	sex
Min. : 2.0	Drug :19	Min. :114.0	Min. :111.0	Min. :100.0	Min. :102.0	Min. : 97.0	Min. :38.00	Length:39
1st Qu.:11.5	Placebo:20	1st Qu.:115.0	1st Qu.:113.0	1st Qu.:112.0	1st Qu.:106.5	1st Qu.:101.5	1st Qu.:42.00	Class :character
Median :21.0	NA	Median :117.0	Median :115.0	Median :113.0	Median :109.0	Median :108.0	Median :48.00	Mode :character
Mean :21.0	NA	Mean :116.7	Mean :114.3	Mean :112.4	Mean :109.4	Mean :106.7	Mean :47.95	NA
3rd Qu.:30.5	NA	3rd Qu.:118.0	3rd Qu.:115.0	3rd Qu.:113.0	3rd Qu.:113.5	3rd Qu.:112.0	3rd Qu.:51.50	NA
Max. :40.0	NA	Max. :121.0	Max. :119.0	Max. :118.0	Max. :117.0	Max. :115.0	Max. :63.00	NA

## [1] "How to process data from wide to long format creating a time variable in one RM variable SCORE is given in the attached code: Long format are desired in eased plotting and mixed models analysis"

About Baseline B0

Baseline B0 value is about the same for the two groups that is SRS Randomization on baseline value is correct, matchingng of covariates seems reasonable (Age ,Sex and range of B0)
Baseline value are slightly higher for older as evidence base physiology (Rigidity artery aging i-e.)
BP increased with aging is proven but in this study is disregarded for study under 4 months time: Hence Age will be a fixed covariates in time.
Baseline value are slightly higher range for women in Placebo.
Conclusion:
As long as baseline is balanced in trt Score value it is fine ! now balancing baseline score within strata is something much more difficult with such low N count (39).

PLOTING LONGITUDINAL DATA

xyplot by id faceted

Both trt result in lowering BP after 4M The magnitude is greater in trt Drug group.
Slope is lower in Placebo group
Linearity in both group might be a crude assumptions
Concave downwards occurred in Drug Group showing an accelerated effect
Inflection point seems to occur at 2 Months
ID 6 & / show erratic/antagonist pattern and might increase variance at that time point.

Spaghetti plot id trajectories

Longitudinal data show often an increased variances in Time due to trt effect, between id difference and within variances: To account for such “type” of variance is not an easy task and terminology between statistician not compatible.

Here we not a first phase a between id variance (“space between point”) B0 and at B2 we notciced some individualization

At time 2-3 the between variance as taken larger values and the 6-7 contribute to the pattern also

At the end of the study we see clustered trt individualization and a stable between patients variance.

Now i will present on lattice same plot faceted by trt:

Spaguetti plot id trajectories | trt

Wiggly pattern occur on both groups

Pronounced slopping is obvious in trt Drug group

Clear group cluster occur after 4 Months. We will try in another blog with MDS techniques to catch that similarities pattern.

Flattening up on the Placebo group is also noticeable that is in polynomial LM we might expect negative vs positive > 1 order coefficients to occurs due to different curvature.

Mean Plot

In the A group a step-down slop occurs after 2 months of active drug Conversely flattening occurs in the Placebo group after 2 months. A linear assumptions might be a crude assumptions but for a accurate modelling ; piece linear of polynomial trend should be envisioned to capture the pattern at 2-3 months periods.

## 'data.frame':    5 obs. of  2 variables:
##  $ Drug   : chr  "116.6842" "113.4211" "110.5789" "106.1053" ...
##  $ Placebo: chr  "116.7500" "115.2000" "114.0500" "112.4500" ...

	Drug	Placebo	difference
baseline	116.6842	116.7500	0
bp1	113.4211	115.2000	-2
bp2	110.5789	114.0500	-4
bp3	106.1053	112.4500	-6
bp4	101.1579	111.9500	-10

1st. Analysis : From Baseline to a 4 Months time point [Naive⁶]

Effect of Drug in lowering is BP at 4 months time is obvious in graphs:
Dotplot + Boxplot show clear clustered drops - no overlap - no outlier
Placebo[B] show no 4th quartile that is Normality assumptions seems doubty
Slope of DRUG A are steeper and more parallel Specila pattern for id 7
Placebo exhibit lower value drop and greater variability in slope :
Note: Keep scaling and range as fixed for both graph Note that on Placebo group several value of Baseline are duplicated despite aspect of that graph check it with command (duplicated(BP$baseline and dp4[BP$trt==“B”]))

Regression Line: comparison at Baseline & only to last time points:

Some figures check to asses some of the assumptions for further inferences:

## [1] "BP baseline value"

##     Drug  Placebo 
## 116.6842 116.7500

## [1] "sd B0 value"

##     Drug  Placebo 
## 2.007297 2.438183

## [1] "Difference from baseline to 4M"

Baseline value mean: around 116.7 both arms: B 116.75 / A 116.68 mmHG:

Note: We don’t know if matching occurred in DOE
SD in baseline is the same that is the one population assumptions

About the research question The effect of TRT?

##     Drug  Placebo 
## 15.56211  4.77000

## [1] "A drop of 15 mmHG might be expected by trt A vs 5 mmHG"

Drop from baseline mean at 4 Month trt mean: B 4.77 / A 15.56 mmHG**

Pooled B0 SD about: 2.2 mmHG

That is:

Drop is greater than resp.pooled sd
High effect size might be expected even in PLA group.
Design experiment plans usually a >= 10 mmHg drop as first clinical intention (i.e Captopril ECA inhibitor).

Inference: T.Test

Although not planned to discuss deeply in this blog (normality etc) we report two type of t.test:

2 sample T test pooled variance at 4Months [_Naive _approach]
DoD T.test pooled variance

Lets graphs first the boxplot of the difference time 4M to baseline value:

The eyes of the “JCVD aware statistician” can tell that if no overlapping( indentation of IQR segment occur a significant difference between arm is likely: Let test it:

## [1] "T-test 4month point"

## 
##  Two Sample t-test
## 
## data:  BP$bp4 by BP$trt
## t = -15.046, df = 37, p-value < 2.2e-16
## alternative hypothesis: true difference in means between group Drug and group Placebo is not equal to 0
## 95 percent confidence interval:
##  -12.245436  -9.338774
## sample estimates:
##    mean in group Drug mean in group Placebo 
##              101.1579              111.9500

Now a clever way to use baseline is to make a DOD T test:

## [1] "T-test DOD difference [Baseline-4months] value score"

## 
##  Two Sample t-test
## 
## data:  BP$baseline - BP$bp4 by BP$trt
## t = 13.203, df = 37, p-value = 1.427e-15
## alternative hypothesis: true difference in means between group Drug and group Placebo is not equal to 0
## 95 percent confidence interval:
##   9.080235 12.372397
## sample estimates:
##    mean in group Drug mean in group Placebo 
##              15.52632               4.80000

## [1] "print making the difference between 2 RV (BBO-B4) remove covariance as pairing  the data within groups:Mistake Often seen: don't put Paired =true between two group as group are independant each other: that is only the RM obesrvations within an ID in group that is paired"

True difference in trt time points 4 two sample T test :

CONFINT: [9.34-12.24]

True difference in DOD T test:

CONFINT: [9.08-12.37]

The interpretation is left at the reader : However note that usually we expect with DOD more precise and difference in the result but here not as much as expected. even less

What could be the explanation?

Looking at correlation / covariance in both arms presented below and which value is below <0.5 might explain why such result occurs.

Differentiating is not removing as much covariance from the DoD test VAR[X-Y] = VARX + VARY - 2 COV XY

This might be a sign also to researcher, that a linear assumptions between baseline and time point 4 might be doubtful due to low correlations (See xyplot by id above) and a CPS (Compound Symmetry structure) of the COV matrix might be doubtful.

Exploring Correlations and variability type:

## [1] "figures now: within groups and all for baseline-4M datapoints"

## [1] NA

## [1] NA

## [1] 0.127816

TRT	COR B0-BP4M
A	0.10
B	0.45
OVERALL	0.12

Note: Usually we expect in each measurements a great variability/instability of BP monthly measure (protocol, devices used)…Here the respective sd for the time variable (each months) is pretty stable and low (less than 5mmHG): That might be a signed that the mean of several BP at each visit were made (?).Researcher has to go back in the design of the experiment to confirm it and have a physiological knowledge of the involved process:

	Drug	Placebo
0	1.945154	2.124419
1	1.677160	1.609184
2	3.271354	1.669384
3	2.558097	2.502104
4	2.007297	2.438183

Linear model

Using BP4 only for modelling [Naive]⁷

Here Baseline will be considered as a Fixed controlled co-variate, not a time co-variate,that is not changing with time.

Note: For sake of short we wont analysis the normality of errors and other assumptions

Note : Note on Contrasts:

To have a possible lowering trend value of the active trt contrast should be done as “TREATMENT” contrast to the reference Placebo group = B as reference level factor.We expect under this kind of contrast a negative regression coefficient for the active drug A if effective.

Remark:

If you want an orthogonal testing lm ANOVA type III, contrast must be set to sum!.

## [1] "confirm type of contrast set and change to PLACEBO 0 as refs"

##         Placebo
## Drug          0
## Placebo       1

##         Drug
## Placebo    0
## Drug       1

ANOVA like TRT vs ANCOVA like TRT +B0⁸

m1 = lm(BP$bp4~BP$trt)

m2 = lm(BP$bp4~BP$trt+BP$baseline)

On simple linear model as TRT is a factor the beta coefficients reflect the decrease in mean in BP from placebo to trt (contr trt ref PLA) and baseline , a continuous covariates represent a slope respectively.

What is special in this presented Technics here :

An Interaction in m2 between baseline(static Xi) and trt as no meaning due to randomization and equal value

- Baseline Wald test and value must not be interpreted.

- Including a baseline as a [usually] a massive impact on SE if cor TRT-B0 > 0.5

## 
## Call:
## lm(formula = BP$bp4 ~ BP$trt)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -4.158 -1.950  0.050  1.446  3.842 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 111.9500     0.5006  223.61   <2e-16 ***
## BP$trtDrug  -10.7921     0.7173  -15.05   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.239 on 37 degrees of freedom
## Multiple R-squared:  0.8595, Adjusted R-squared:  0.8557 
## F-statistic: 226.4 on 1 and 37 DF,  p-value: < 2.2e-16

## 
## Call:
## lm(formula = BP$bp4 ~ BP$trt + BP$baseline)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.9210 -1.6971 -0.2505  1.5233  4.0676 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  73.4793    20.3855   3.604 0.000939 ***
## BP$trtDrug  -10.7704     0.6937 -15.525  < 2e-16 ***
## BP$baseline   0.3295     0.1746   1.888 0.067149 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.165 on 36 degrees of freedom
## Multiple R-squared:  0.8722, Adjusted R-squared:  0.8651 
## F-statistic: 122.8 on 2 and 36 DF,  p-value: < 2.2e-16

\[ \operatorname{\widehat{BP\$bp4}} = 111.95 - 10.79() \] \[ \operatorname{\widehat{BP\$bp4}} = 73.48 - 10.77() + 0.33(\operatorname{BP\$baseline}) \]

## [1] "confidence intervals of coeeficients m1 vs m2"

##              2.5 % 97.5 %
## (Intercept) 110.94 112.96
## BP$trtDrug  -12.25  -9.34

##              2.5 % 97.5 %
## (Intercept)  32.14 114.82
## BP$trtDrug  -12.18  -9.36
## BP$baseline  -0.02   0.68

Fitted value

## [1] "plotm1 fitted values"

## [1] "plotm2 fitted values"

TRT is highly significant in both models.

Recall:On simple linear model with TRT is a factor the beta coefficients reflect the decrease in mean in BP from placebo (contrast trt ref PLA)

The SE of the TRT including a baseline is a of little gain in precision hence Inference but usually when higher correlation TRT -B= exist SE are far more impressive.

Note the T value for baseline is the same for anova model comparison

Note the Anova model comparison is a limited value here anova(m1,m2)

A little gain is done here in this dataset with but with prediction points of view a lot is done by adding a continuous covariates (ANCOVA sensus.lato) vs an ANOVA: lm as E(Y|X) is always better that the mean.

However the prediction of the drug A is fairly in accordance with true value.

Anova model comparison TYPE I vs III unbalanced design

We notice a slight but negligible difference in type I and III due to missing ID1 unbalanced:

## Anova Table (Type III tests)
## 
## Response: BP$bp4
##              Sum Sq Df  F value    Pr(>F)    
## (Intercept)   60.91  1  12.9924  0.000939 ***
## BP$baseline   16.71  1   3.5634  0.067149 .  
## BP$trt      1129.97  1 241.0299 < 2.2e-16 ***
## Residuals    168.77 36                       
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

##             Df Sum Sq Mean Sq F value Pr(>F)    
## BP$trt       1 1134.8  1134.8 242.067 <2e-16 ***
## BP$baseline  1   16.7    16.7   3.563 0.0671 .  
## Residuals   36  168.8     4.7                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

LM with Y = Difference from baseline = Y~BP0

Now lets explore an DOD Y | X = (BP4- BPbaseline) “like” regression with & without baseline in the X covariate

m3 = lm(BP$bp4-BP$Baseline~BP$trt)

m2 = lm(BP$bp4-BP$baseline~BP$trt+BP$baseline)

## 
## Call:
## lm(formula = BP$diff4b ~ BP$trt)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -5.474 -1.837 -0.200  1.526  5.800 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  -4.8000     0.5670  -8.465 3.50e-10 ***
## BP$trtDrug  -10.7263     0.8124 -13.203 1.43e-15 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.536 on 37 degrees of freedom
## Multiple R-squared:  0.8249, Adjusted R-squared:  0.8202 
## F-statistic: 174.3 on 1 and 37 DF,  p-value: 1.427e-15

## 
## Call:
## lm(formula = BP$diff4b ~ BP$baseline + BP$trt)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.9210 -1.6971 -0.2505  1.5233  4.0676 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  73.4793    20.3855   3.604 0.000939 ***
## BP$baseline  -0.6705     0.1746  -3.841 0.000478 ***
## BP$trtDrug  -10.7704     0.6937 -15.525  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.165 on 36 degrees of freedom
## Multiple R-squared:  0.8758, Adjusted R-squared:  0.8689 
## F-statistic: 126.9 on 2 and 36 DF,  p-value: < 2.2e-16

##              2.5 % 97.5 %
## (Intercept)  -5.95  -3.65
## BP$trtDrug  -12.37  -9.08

##              2.5 % 97.5 %
## (Intercept)  32.14 114.82
## BP$baseline  -1.02  -0.32
## BP$trtDrug  -12.18  -9.36

\[ \operatorname{\widehat{BP\$bp4}} = 111.95 - 10.79() \] \[ \operatorname{\widehat{BP\$bp4}} = 73.48 - 10.77() + 0.33(\operatorname{BP\$baseline}) \]

Fitted values

Here nothing is really gain in R2 squared including a baseline in Xis (but nothing is really lost too)!

No co-linearity is observed in the hat matrix because it depends only on X (VIF ::car Variance inflation factor): However we usually check co-linearity with 2 continuous variables but including baseline left and right from the equation is not ill conditioned matrix result.

## [1] "VIF m4"

## BP$baseline      BP$trt 
##    1.000274    1.000274

Note : At this stage the proposal of robust estimator is not of interest no repeated measures designed are made here on Yis even subtraction baseline (correlation is removed).

Which simple lm model should we choose m2 vs m4 ? :

Lets summarise the model performance

## 
## Call:
## lm(formula = BP$bp4 ~ BP$trt + BP$baseline)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.9210 -1.6971 -0.2505  1.5233  4.0676 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  73.4793    20.3855   3.604 0.000939 ***
## BP$trtDrug  -10.7704     0.6937 -15.525  < 2e-16 ***
## BP$baseline   0.3295     0.1746   1.888 0.067149 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.165 on 36 degrees of freedom
## Multiple R-squared:  0.8722, Adjusted R-squared:  0.8651 
## F-statistic: 122.8 on 2 and 36 DF,  p-value: < 2.2e-16

## 
## Call:
## lm(formula = BP$diff4b ~ BP$baseline + BP$trt)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.9210 -1.6971 -0.2505  1.5233  4.0676 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  73.4793    20.3855   3.604 0.000939 ***
## BP$baseline  -0.6705     0.1746  -3.841 0.000478 ***
## BP$trtDrug  -10.7704     0.6937 -15.525  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.165 on 36 degrees of freedom
## Multiple R-squared:  0.8758, Adjusted R-squared:  0.8689 
## F-statistic: 126.9 on 2 and 36 DF,  p-value: < 2.2e-16

Well R2 sq are similar , RMSE the same:

As advice choose the m2 model lm BP4~Baseline + trt in the original scale that the one mostly used in preliminary analysis.[see Zhang,2014]

Preliminary conclusions :

Including baseline as covariate (X) improved definitely the prediction of your model (R2sq, SE of Coeff) especially if correlation is above 0.5.

LONGITUDINAL DATA ANALYSIS

Using the full information of time

Let first explore the possibility of an interactions time * trt effect: As mentioned Interactions , that is slope in time with trt can be visualized as such on the interaction plot. Time*trt effect show how the magnitude of drugs acts in time (and induce a increase in overall variance )

Therefore it might be assessed that trt difference between arms is effective in time elapsed: (fan out sloping)

Before you dataset must be revamped in long format and adding a static baseline in long format

id	trt	age	sex	bp	Score
2	Drug	51	M	baseline	116
3	Drug	48	F	baseline	119
4	Drug	42	F	baseline	115
5	Drug	49	M	baseline	116
6	Drug	47	M	baseline	117
7	Drug	50	F	baseline	118

Lets explore a basic longitudinal model where baseline is included as a dynamic time points:

(ml1 Y~only time)
ml2 Y~time+trt
ml3 time*trt (ANCOVA)

## [1] "confirm type of contrast set and change to PLACEBO 0 as refs"

##         Drug
## Placebo    0
## Drug       1

## [1] "set PLACEBO=B AS REF [0;0]"

##         Drug
## Placebo    0
## Drug       1

## [1] "ml1"

## 
## Call:
## lm(formula = Score ~ time, data = new)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -11.8923  -2.3897   0.1077   2.1051   8.1128 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 116.8974     0.4560  256.37   <2e-16 ***
## time         -2.5026     0.1861  -13.44   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.676 on 193 degrees of freedom
## Multiple R-squared:  0.4836, Adjusted R-squared:  0.4809 
## F-statistic: 180.7 on 1 and 193 DF,  p-value: < 2.2e-16

## [1] "ml2"

## 
## Call:
## lm(formula = Score ~ time + trt, data = new)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -9.5895 -1.8347 -0.0851  2.4054  5.9251 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 119.0851     0.4142  287.53   <2e-16 ***
## time         -2.5026     0.1474  -16.98   <2e-16 ***
## trtDrug      -4.4905     0.4169  -10.77   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.91 on 192 degrees of freedom
## Multiple R-squared:  0.6781, Adjusted R-squared:  0.6747 
## F-statistic: 202.2 on 2 and 192 DF,  p-value: < 2.2e-16

## [1] "ml3"

## 
## Call:
## lm(formula = Score ~ time * trt, data = new)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -9.5895 -1.4263 -0.2632  1.5737  4.4500 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  116.5500     0.3896 299.174  < 2e-16 ***
## time          -1.2350     0.1590  -7.765 4.84e-13 ***
## trtDrug        0.7132     0.5581   1.278    0.203    
## time:trtDrug  -2.6018     0.2279 -11.419  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.249 on 191 degrees of freedom
## Multiple R-squared:  0.8087, Adjusted R-squared:  0.8057 
## F-statistic: 269.1 on 3 and 191 DF,  p-value: < 2.2e-16

## [1] "conf:int ml1,2,3"

##                  2.5 %     97.5 %
## (Intercept) 115.998108 117.796764
## time         -2.869713  -2.135415

##                  2.5 %     97.5 %
## (Intercept) 118.268236 119.902020
## time         -2.793201  -2.211927
## trtDrug      -5.312843  -3.668210

##                    2.5 %      97.5 %
## (Intercept)  115.7815837 117.3184163
## time          -1.5487046  -0.9212954
## trtDrug       -0.3877528   1.8140685
## time:trtDrug  -3.0512870  -2.1523972

## [1] "anova ssq ml1,2,3"

## Analysis of Variance Table
## 
## Model 1: Score ~ time
## Model 2: Score ~ time + trt
## Model 3: Score ~ time * trt
##   Res.Df     RSS Df Sum of Sq      F    Pr(>F)    
## 1    193 2608.24                                  
## 2    192 1625.85  1    982.39 194.19 < 2.2e-16 ***
## 3    191  966.25  1    659.60 130.38 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

## [1] "sigma ml1,2,3"

## [1] 3.676163

## [1] 2.909975

## [1] 2.249196

ml1 as no meaning in RCT as TRT in not tested

We definitely see that ml3 is the convenient model (sigma,Rsquared) and is in accordance with the plot we saw preceedingly: In ANCOVA (contarst trt) the negative sign indicated that drugs reduce further each month the value of Drug patient of 2.60 mmHG in (linear) mean.

Note In ANCOVA trt factor must be conserved in the OLS but is not of interpretation (see Fox,2015).

See seq ANOVA for factor trt < 0.05 in ml3 :

## Analysis of Variance Table
## 
## Response: Score
##            Df  Sum Sq Mean Sq F value    Pr(>F)    
## time        1 2442.50 2442.50  482.81 < 2.2e-16 ***
## trt         1  982.39  982.39  194.19 < 2.2e-16 ***
## time:trt    1  659.60  659.60  130.38 < 2.2e-16 ***
## Residuals 191  966.25    5.06                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

The model matrix of ml3 look like this:

(Intercept)	time	trtDrug	time:trtDrug
1	0	1	0
1	0	1	0
1	0	1	0
1	0	1	0
1	0	1	0
1	0	1	0
1	0	1	0
1	0	1	0

Moreover we observed in the interaction plot on non equal interval spacing in time but also on the spaguetti plot by ID.

One way to test this is by changing the time as an ordered factor: “factor(time,ordered=TRUE”:

## 
## Call:
## lm(formula = Score ~ factor(time, ordered = TRUE) * trt, data = new)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -10.5789  -1.4211  -0.1579   1.5789   4.5500 
## 
## Coefficients:
##                                         Estimate Std. Error t value Pr(>|t|)
## (Intercept)                            114.08000    0.22318 511.154  < 2e-16
## factor(time, ordered = TRUE).L          -3.90541    0.49905  -7.826 3.76e-13
## factor(time, ordered = TRUE).Q           0.44098    0.49905   0.884   0.3780
## factor(time, ordered = TRUE).C           0.22136    0.49905   0.444   0.6579
## factor(time, ordered = TRUE)^4           0.28685    0.49905   0.575   0.5661
## trtDrug                                 -4.49053    0.31975 -14.044  < 2e-16
## factor(time, ordered = TRUE).L:trtDrug  -8.22775    0.71499 -11.508  < 2e-16
## factor(time, ordered = TRUE).Q:trtDrug  -1.77729    0.71499  -2.486   0.0138
## factor(time, ordered = TRUE).C:trtDrug  -0.50430    0.71499  -0.705   0.4815
## factor(time, ordered = TRUE)^4:trtDrug   0.09688    0.71499   0.135   0.8924
##                                           
## (Intercept)                            ***
## factor(time, ordered = TRUE).L         ***
## factor(time, ordered = TRUE).Q            
## factor(time, ordered = TRUE).C            
## factor(time, ordered = TRUE)^4            
## trtDrug                                ***
## factor(time, ordered = TRUE).L:trtDrug ***
## factor(time, ordered = TRUE).Q:trtDrug *  
## factor(time, ordered = TRUE).C:trtDrug    
## factor(time, ordered = TRUE)^4:trtDrug    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.232 on 185 degrees of freedom
## Multiple R-squared:  0.8176, Adjusted R-squared:  0.8087 
## F-statistic: 92.11 on 9 and 185 DF,  p-value: < 2.2e-16

Growth model

To asses what type of polynomial you need depends also on your goal but as seen on graph some concave up and concave down pattern might even required a third degree polynomial?!

## [1] "model polynomial 3 order in time"

## 
## Call:
## lm(formula = Score ~ poly(time, 3), data = new)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -12.119  -2.493   0.322   2.101   8.348 
## 
## Coefficients:
##                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)    111.8923     0.2643 423.399   <2e-16 ***
## poly(time, 3)1 -49.4217     3.6904 -13.392   <2e-16 ***
## poly(time, 3)2  -2.6534     3.6904  -0.719    0.473    
## poly(time, 3)3  -0.1519     3.6904  -0.041    0.967    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.69 on 191 degrees of freedom
## Multiple R-squared:  0.485,  Adjusted R-squared:  0.4769 
## F-statistic: 59.96 on 3 and 191 DF,  p-value: < 2.2e-16

## [1] "model polynomial order 2 in the interactions"

## 
## Call:
## lm(formula = Score ~ poly(time, 2) * trt, data = new)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -10.3038  -1.5489  -0.1971   1.2729   4.2729 
## 
## Coefficients:
##                        Estimate Std. Error t value Pr(>|t|)    
## (Intercept)            114.0800     0.2216 514.714  < 2e-16 ***
## poly(time, 2)1         -24.3893     3.0950  -7.880 2.51e-13 ***
## poly(time, 2)2           2.7539     3.0950   0.890   0.3747    
## trtDrug                 -4.4905     0.3175 -14.142  < 2e-16 ***
## poly(time, 2)1:trtDrug -51.3823     4.4342 -11.588  < 2e-16 ***
## poly(time, 2)2:trtDrug -11.0992     4.4342  -2.503   0.0132 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.216 on 189 degrees of freedom
## Multiple R-squared:  0.8162, Adjusted R-squared:  0.8113 
## F-statistic: 167.8 on 5 and 189 DF,  p-value: < 2.2e-16

## [1] "anova"

## Analysis of Variance Table
## 
## Model 1: Score ~ poly(time, 2) * trt
## Model 2: Score ~ poly(time, 3)
##   Res.Df     RSS Df Sum of Sq      F    Pr(>F)    
## 1    189  928.43                                  
## 2    191 2601.17 -2   -1672.7 170.26 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Probably a reasonable assumptions for other linear model is an order 2 polynomial with trt poly order 2 interaction to account for change in simple curvature change up/down.

Until then we haven’t accounted for the RM correlations that is the SE are inflated and inference not precised and should not be used for publications.

For those who are not used to LMM you can use a gls function and ask for a Robust Sandwich :: robust :: estimator of regression coefficients for the retained model (See Fox,2015),

Including B0 as static covariate

Now will explore that baseline B0 is a static predictor not included in the dynamic time frame Bp1,Bp2,Bp3,Bp4. Therefore for each line of the model matix the baseline B0 for each id must be repeated: This required a little bit involvement in reshaping the data in R

## [1] "X matrix and dimension of restructed data"

##    id  trt age sex baseline  bp Score
## 1   2 Drug  51   M      116 bp1   113
## 2   3 Drug  48   F      119 bp1   115
## 3   4 Drug  42   F      115 bp1   113
## 4   5 Drug  49   M      116 bp1   112
## 5   6 Drug  47   M      117 bp1   112
## 6   7 Drug  50   F      118 bp1   111
## 7   8 Drug  61   M      120 bp1   115
## 8   9 Drug  43   M      114 bp1   112
## 9  10 Drug  51   M      115 bp1   113
## 10 11 Drug  47   F      117 bp1   112

## [1] 156   7

## [1] "X matrix and dimension of long format bp0,1,2,3,4 preceeding mdel datadata "

##    id  trt age sex       bp Score time
## 1   2 Drug  51   M baseline   116    0
## 2   3 Drug  48   F baseline   119    0
## 3   4 Drug  42   F baseline   115    0
## 4   5 Drug  49   M baseline   116    0
## 5   6 Drug  47   M baseline   117    0
## 6   7 Drug  50   F baseline   118    0
## 7   8 Drug  61   M baseline   120    0
## 8   9 Drug  43   M baseline   114    0
## 9  10 Drug  51   M baseline   115    0
## 10 11 Drug  47   F baseline   117    0

## [1] 195   7

## [1] "BP must be a time numeric value"

##    id  trt age sex baseline  bp Score time
## 1   2 Drug  51   M      116 bp1   113    1
## 2   3 Drug  48   F      119 bp1   115    1
## 3   4 Drug  42   F      115 bp1   113    1
## 4   5 Drug  49   M      116 bp1   112    1
## 5   6 Drug  47   M      117 bp1   112    1
## 6   7 Drug  50   F      118 bp1   111    1
## 7   8 Drug  61   M      120 bp1   115    1
## 8   9 Drug  43   M      114 bp1   112    1
## 9  10 Drug  51   M      115 bp1   113    1
## 10 11 Drug  47   F      117 bp1   112    1

## [1] "confirm type of contrast set and change to PLACEBO 0 as refs"

##         Drug
## Placebo    0
## Drug       1

## [1] "set PLACEBO=B AS REF"

## [1] "note If you messed up in cyour relvel factor in R interetation could be catastrophically wrong"

##         Drug
## Placebo    0
## Drug       1

## [1] "Lets report the coefficient coefficients confint"

##                   2.5 %     97.5 %
## (Intercept)  58.7708090 99.8081592
## time         -1.5702558 -0.6997442
## trtDrug       0.1945924  3.6102205
## baseline      0.1411263  0.4920303
## time:trtDrug -3.6149071 -2.3677245

The retained model with a baseline as static predictor tell us that every month a linear difference of -3 mmHG for the drug occurs that is for a 4 Months treatment will be about 12 mmHG *[ 14.4 : 9.48* ] of E(Y|X) difference will occur:Lets check it with table of difference form the mean difference at time 4:

	Drug	Placebo	difference
baseline	116.6842	116.7500	0
bp1	113.4211	115.2000	-2
bp2	110.5789	114.0500	-4
bp3	106.1053	112.4500	-6
bp4	101.1579	111.9500	-10

And from the the DoD TTest

9.08 12.37

The linear model gives higher difference result and for a DRUG that is the aim for!

You have to recall that in linear model the E(Y|X) is not the mean of group: It integrate the conditional linear trend within all months of experiment that the real power of regression MODEL.

In Anova the empirical group mean is used instead.

Of course now the main effect of time and trt in our lm cannot be interpreted promptly as J.Fox reported in 2015.

But in bio-statistics what is crucial is how much difference between group is accounted in time (each month) : Therefore the interaction term is the major focus on the model.

Now let see which of our two lm with B0 as time predictor and the other as static predictor is different:

## [1] "model with B0 as dynamic predictor"

## 
## Call:
## lm(formula = Score ~ time * trt, data = new)
## 
## Coefficients:
##  (Intercept)          time       trtDrug  time:trtDrug  
##     116.5500       -1.2350        0.7132       -2.6018

## [1] "effect is 10.40 MMHG drop in favour of DRUG"

## [1] "model with B0 as controlled X static predictor"

## 
## Call:
## lm(formula = Score ~ time * trt + baseline, data = newb)
## 
## Coefficients:
##  (Intercept)          time       trtDrug      baseline  time:trtDrug  
##      79.2895       -1.1350        1.9024        0.3166       -2.9913

## [1] "confint"

##               2.5 % 97.5 %
## (Intercept)  115.78 117.32
## time          -1.55  -0.92
## trtDrug       -0.39   1.81
## time:trtDrug  -3.05  -2.15

##              2.5 % 97.5 %
## (Intercept)  58.77  99.81
## time         -1.57  -0.70
## trtDrug       0.19   3.61
## baseline      0.14   0.49
## time:trtDrug -3.61  -2.37

In the first model ml3 : Score~time*trt the interaction effect is 10.40 mmHG 4*[3.05:2.15MMHG drop in favor of DRUG after 4 months-

In the first model ml4B0 the interaction effect is about 12 mmHG 4*[3.61;2.37] drop in favor of DRUG after 4 months-

Conclusions

Including a baseline as a static controlled covariate X result in better efficiency in effect and more precise assessment of a treatment confidence intervall in RCT

The wrong way to do it

Wrong way to do it: Wide and bp as factor time!

## 
## Call:
## lm(formula = BP$bp4 ~ BP$bp1 + BP$bp3 + BP$baseline + BP$trt)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.7508 -1.3071 -0.3603  1.4835  5.0010 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  28.1747    25.4997   1.105    0.277    
## BP$bp1        0.2199     0.2498   0.881    0.385    
## BP$bp3        0.3484     0.1328   2.623    0.013 *  
## BP$baseline   0.1650     0.1971   0.837    0.409    
## BP$trtDrug   -8.1794     1.0747  -7.611 7.64e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.983 on 34 degrees of freedom
## Multiple R-squared:  0.8988, Adjusted R-squared:  0.8869 
## F-statistic: 75.47 on 4 and 34 DF,  p-value: < 2.2e-16

Linear Mixed Model will be presented in our next BLOG

STAY tuned

REFERENCES

Chen, D. (., Peace, K. E. (2011). Clinical Trial Data Analysis Using R. États-Unis: CRC Press.

Fox, J., & Weisberg, S. (2018). An R companion to applied regression. Sage publications.

Zhang, S., Paul, J.E., Nantha-Aree, M., Buckley, N., Shahzad, U., Cheng, J., Debeer, J., Winemaker, M.J., Wismer, D., Punthakee, D., Avram, V., & Thabane, L. (2014). Empirical comparison of four baseline covariate adjustment methods in analysis of continuous outcomes in randomized controlled trials. Clinical Epidemiology, 6, 227 - 235.

A first approcach would be to use a robust estimateor of the COV matrix knwon as the sandwich / Hubert estimators↩︎
- or Growth curve if ciurvature pattern
↩︎
Large clinical trials have shown that combined treatment with reserpine plus a thiazide diuretic reduces mortality of people with hypertension..↩︎
Back on tracking the history of this RCT (part of the job) it is not clear from source which treatment was tested either Resepinic or betaB…The intention here is to demonstrate the technical statistical aspect with repeated data not a biophysical full proof of a drug per se.↩︎
Missing value will be replace by a Propensity score matching approach in this slight unbalanced design to make an repeated Anova and result compared.↩︎
Naive because it doesn’t take into account repeated corrrelated measures hence partial information only Use LMM or GEE↩︎
Because we use only one type of information at time 4! some study are so designed but be aware that slope information might be misleading if one to compare the baseline with each trt with low sample size such exposed in Captopril se xyplot study↩︎
[Naive model] doesn’t take accout of RM.↩︎

Longitudinal regression and Baseline inclusion

Mudryjm

2025-11-21

GOAL

RCT and DATA

Know your Data: Systematic review

Data summary

About Baseline B0

PLOTING LONGITUDINAL DATA

xyplot by id faceted

Spaghetti plot id trajectories

Spaguetti plot id trajectories | trt

Mean Plot

1st. Analysis : From Baseline to a 4 Months time point [Naive⁶]

Regression Line: comparison at Baseline & only to last time points:

About the research question The effect of TRT?

Inference: T.Test

Exploring Correlations and variability type:

Linear model

Using BP4 only for modelling [Naive]⁷

ANOVA like TRT vs ANCOVA like TRT +B0⁸

LM with Y = Difference from baseline = Y~BP0

Which simple lm model should we choose m2 vs m4 ? :

Preliminary conclusions :

LONGITUDINAL DATA ANALYSIS

Using the full information of time

Growth model

Including B0 as static covariate

Conclusions

The wrong way to do it

REFERENCES

Longitudinal regression and Baseline inclusion

Mudryjm

2025-11-21

GOAL

RCT and DATA

Know your Data: Systematic review

Data summary

About Baseline B0

PLOTING LONGITUDINAL DATA

xyplot by id faceted

Spaghetti plot id trajectories

Spaguetti plot id trajectories | trt

Mean Plot

1st. Analysis : From Baseline to a 4 Months time point *[Naive*6]

Regression Line: comparison at Baseline & only to last time points:

About the research question The effect of TRT?

Inference: T.Test

Exploring Correlations and variability type:

Linear model

Using BP4 only for modelling [Naive]7

ANOVA like TRT vs ANCOVA like TRT +B08

LM with Y = Difference from baseline = Y~BP0

Which simple lm model should we choose m2 vs m4 ? :

Preliminary conclusions :

LONGITUDINAL DATA ANALYSIS

Using the full information of time

Growth model

Including B0 as static covariate

Conclusions

The wrong way to do it

REFERENCES

1st. Analysis : From Baseline to a 4 Months time point [Naive⁶]

Using BP4 only for modelling [Naive]⁷

ANOVA like TRT vs ANCOVA like TRT +B0⁸