1 Description

The data come from a randomized, double-blind clinical trial comparing active hypnotic drug with a placebo in patients who have insomnia problems. The response is the patient’s reported time (in minutes) to fall asleep after going to bed. Patients responded before and after a two-week treatment period. The patients receiving the two treatments were independent samples.

Column 1: Subject ID
Column 2: Treatment, Active = 1, Placebo = 0
Column 3: Time, Initial = 0, Follow-up = 1
Column 4: Time to sleep, < 20 min. = 1, 20-30 = 2, 30-60 = 3, > 60 min. = 4

2 input data

dta2 <- read.table("C:/Users/HANK/Desktop/HOMEWORK/insomnia.txt", h=T)
# rename
colnames(dta2) <- c("ID", "Treatment", "Time", "Sleeplatency")
head(dta2)
##   ID Treatment Time Sleeplatency
## 1  1         1    0            1
## 2  1         1    1            1
## 3  2         1    0            1
## 4  2         1    1            1
## 5  3         1    0            1
## 6  3         1    1            1
str(dta2)
## 'data.frame':    478 obs. of  4 variables:
##  $ ID          : int  1 1 2 2 3 3 4 4 5 5 ...
##  $ Treatment   : int  1 1 1 1 1 1 1 1 1 1 ...
##  $ Time        : int  0 1 0 1 0 1 0 1 0 1 ...
##  $ Sleeplatency: int  1 1 1 1 1 1 1 1 1 1 ...

3 GEE

# repolr [data = as.data.frame()]
pacman::p_load(repolr)
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
# repolr
summary(dta2_repolr <- repolr(Sleeplatency ~ Treatment + Time + Treatment:Time,
                 subjects = "ID",
                 data = as.data.frame(dta2),
                 times = c(1,2),
                 categories = 4,
                 corr.mod = "independence",
                 po.test=TRUE,))
## 
## repolr: 2016-02-26 version 3.4 
## 
## Call:
## repolr(formula = Sleeplatency ~ Treatment + Time + Treatment:Time, 
##     subjects = "ID", data = as.data.frame(dta2), times = c(1, 
##         2), categories = 4, corr.mod = "independence", po.test = TRUE)
## 
## Coefficients: 
##                 coeff     se.robust  z.robust  p.value 
## cuts1|2          -2.2671    0.2091   -10.8422    0.0000
## cuts2|3          -0.9515    0.1769    -5.3787    0.0000
## cuts3|4           0.3517    0.1751     2.0086    0.0446
## Treatment         0.0336    0.2369     0.1418    0.8872
## Time              1.0381    0.1564     6.6375    0.0000
## Treatment:Time    0.7077    0.2458     2.8792    0.0040
## 
## Correlation Structure:  independence 
## Fixed Correlation:  0 
## PO Score Test: 8.6394 (d.f. = 6 and p.value = 0.1949)
# polr

sjPlot::tab_model(m0 <- MASS::polr(ordered(Sleeplatency) ~ Treatment + Time + 
                                   Treatment:Time, data=dta2, Hess = T),
                  show.obs=F, show.r2=F, show.se=T, transform=NULL)
  ordered(Sleeplatency)
Predictors Log-Odds std. Error CI p
1|2 -2.27 0.20 -2.73 – -1.80 <0.001
2|3 -0.95 0.18 -1.43 – -0.48 <0.001
3|4 0.35 0.17 -0.30 – 1.01 0.045
Treatment -0.03 0.24 -0.50 – 0.43 0.888
Time -1.04 0.24 -1.51 – -0.57 <0.001
Treatment * Time -0.71 0.33 -1.36 – -0.05 0.035 
summary(dta2_repolr$gee)
## Length  Class   Mode 
##      0   NULL   NULL

4 hand-made cumulative prop. against time to sleep

dta2 <- dta2 %>% 
  mutate(Treatment = factor(Treatment, 
                            levels=0:1, 
                            labels=c("Placebo","Active")),
         Time = factor(Time, 
                       levels=0:1, 
                       labels=c("Initial","Follow-up")),
         Sleeplatency = factor(Sleeplatency))

dta2.tbl <- ftable(dta2$Treatment, dta2$Time, dta2$Sleeplatency)

dta2.tbl
##                     1  2  3  4
##                               
## Placebo Initial    14 20 35 51
##         Follow-up  31 29 35 25
## Active  Initial    12 20 40 47
##         Follow-up  40 49 19 11

兩組的睡眠時間均得到改善。其中，採用催眠藥物組(active hypnotic drug)比安慰劑組效果更好。

5 Cumulative prop.plot

obsp <- as.numeric(t(prop.table(dta2.tbl, m=1)))

clp <- c(cumsum(obsp[1:4]), 
         cumsum(obsp[5:8]),
         cumsum(obsp[9:12]), 
         cumsum(obsp[13:16]))

yt <- rep(c(20,25,45,60),4)

y <- data.frame(time2sleep=yt, clprop=clp)

plot(y[,1], y[,2], type="n",
     xlab="Time to Sleep (min.)", ylab="Cumulative proportion")

lines(y[1:4,  1],  y[1:4,2], lty=1, col="steelblue")
lines(y[5:8,  1],  y[5:8,2], lty=1, col="steelblue")
lines(y[9:12, 1],  y[9:12,2], lty=2, col="indianred")
lines(y[13:16,1],  y[13:16,2], lty=2, col="indianred")

text(22, .12, "T,I")
text(22, .47, "T,F")
text(22, .20, "P,I")
text(22, .33, "P,F")

grid()

採用催眠藥物組和安慰劑組間的比較結果:起初皆以相似的時間入睡。但在兩週後，接受催眠藥物治療的人比起安慰劑組，能更迅速入眠。

6 The end