-
Notifications
You must be signed in to change notification settings - Fork 0
/
report.Rnw
262 lines (245 loc) · 19 KB
/
report.Rnw
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
\documentclass[11pt,a4paper]{article}
\title{ASME assignment}
\author{Student number: 106936}
\usepackage{booktabs}
\usepackage{xspace}
\usepackage{amsmath}%
\usepackage{amsfonts}
\usepackage{amssymb}
%\usepackage{graphicx}
\usepackage[bottom=2cm,margin=2.5cm]{geometry}
\usepackage{lmodern}
\usepackage{rotating}
\usepackage[T1]{fontenc}
\usepackage{fancyhdr}
\fancyhead{}
\fancyfoot{}
\setlength{\headheight}{15.2pt}
%\setlength{\footskip=20pt}
\pagestyle{fancy}
\lhead[\thepage]{Student 106936}
\chead[ ASME assignment]{ASME assignment}
\rhead[ Student 106936]{\thepage}
\usepackage[compact]{titlesec}
\titleformat*{\section}{\bfseries}
\titlespacing{\section}{0pt}{1ex}{0ex}
\makeatletter
\newcommand\gobblepars{%
\@ifnextchar\par%
{\expandafter\gobblepars\@gobble}%
{}}
\makeatother
\begin{document}
%\SweaveOpts{concordance=TRUE}
<<setup, results='asis', echo=FALSE, message=FALSE, warning=FALSE>>=
library(survival)
library(foreign)
library(plyr)
library(xtable)
# functions ####
perc <- function(d,n){
round((d/n)*100,2)
}
rate_per_k <- function(x,y){
(x/y)*1000
}
srOrWrapper <- function(x){
cbind(HR=round(1/exp(coef(x)), 2),
se=round(summary(x)$table[,2], 2),
lci=round((1/exp(coef(x)))/exp(1.96*summary(x)$table[,2]), 2),
uci=round((1/exp(coef(x)))*exp(1.96*summary(x)$table[,2]), 2),
p=round(summary(x)$table[,4], 4))
}
se_log_rr <- function(x,y){
sqrt((1/x)+(1/y))
}
z_rr <- function(x,y){
log(x)/y
}
p_rr <- function(x){
2*pnorm(-abs(x))
}
decimal2<-function(x){paste(x,".00",sep="")} # decimal replacer for character strings, 2 decimal places
decimal1<-function(x){paste(x,"0",sep="")} # decimal replacer for character strings, 1 decimal place
simpleCap<-function(x){
s<-as.character(x)
s<-paste(toupper(substring(s,1,1)),substring(s,2),sep="")
# s<-as.factor(s)
}
extractCI <- function(m){
beta <- round(coef(m)["coef"], 3)
se <- round(coef(m)["se(coef)"], 3)
p <- round(coef(m)["Pr(>|z|)"], 3)
ci <- round(m$conf.int[c("lower .95", "upper .95")],3)
res <- rbind(beta, round(exp(beta),3), ci[1], ci[2], p)
res
}
# reading data in ####
first <- read.dta("E:\\My Documents B\\MSc stuff\\SME\\asme\\firstrv.dta")
#first<-read.dta("/home/simon/Documents/MSc_modules/asme/ASMEdata/firstrv.dta")
#first<-read.dta("/home/simon/Documents/MSc/Modules/ASME/firstrv.dta")
# can start with very basics. How many children, how many events. Age and sex.
first$lbw <- as.character(first$lbw)
all <- read.dta("E:\\My Documents B\\MSc stuff\\SME\\asme\\allrv.dta")
#all <- read.dta("/home/simon/Documents/MSc_modules/asme/ASMEdata/allrv.dta")
#all <- read.dta("/home/simon/Documents/MSc/Modules/ASME//allrv.dta")
# get max number of episodes by id and time at risk
all$case <- 0
all$case[all$exittype == "RV diarrhoea"] <- 1
all$dar <- as.numeric(all$endfoll - all$startfoll,format="days")
all$doe2 <- as.numeric(all$startfoll, format="days")
all$dexit <- as.numeric(all$endfoll, format="days")
max.episodes <- ddply(all, .(id), summarise, max=max(numprevepisode, na.rm=TRUE), tdar=sum(dar, na.rm=TRUE))
# Combine max episodes and tdar into first dataframe
first <- merge(first, max.episodes, all.x=TRUE, by.x="id", by.y="id")
first$tyar <- first$tdar/365.25
first$doe2 <- as.numeric(first$doe)
first$dexit <-as.numeric(first$exitdate)
first$fu.time <- first$dexit - first$doe2
first$case <- 0
first$case[first$exittype == "RV diarrhoea"] <- 1
# end read in ####
@
\section{Methods}
Associations between potential risk factors and the rate of rotaviral diarrhoea was examined first by calculating rate ratios and tested using \textit{z}-tests. \\
\indent Interactions between variables were examined first by Mantel-Haenszel analyses where interactions were considered important if they had a substantial effect on the rate ratio.
Important interactions in such analyses were then examined by inclusion or exclusion in multivariable models and were considered important if they had a substantial effect on the rate ratio.\\
\indent To adjust for the effect of calendar time and age, data were split on these variables.
Since the first few months of age were considered likely to show most variation in the rate of rotaviral diarrhoea, age was split into bands of increasing width.
Preliminary analyses suggested that the highest rates of rotaviral diarrhoea occurred between December and February and in June and July.
Calendar time was therefore split by quarter, starting from December.
The period with the greatest person-time was chosen as the baseline period.\\
\indent All variables were considered potential risk factors \textit{a priori}.
Both being of a low birth weight and history of neonatal rotaviral diarrhoea could predispose to rotaviral disease later.
Mother's education level and a family's socio-economic status (SES) could reflect a number of aspects about the environment in which subjects lived and a lower SES was presumed to increase the risk of rotaviral diarrhoea.
Household size would be likely to increase risk of rotaviral diarrhoea by increasing the number of close contacts a child had.
Animal ownership was expected to increase the risk of disease by decreasing general hygiene.
Beediwork could also decrease the general level of hygiene within a household.
Therefore, maximal regression models were fitted.
Poisson regression was chosen over Cox regression as it offers direct estimation of a rate.
To account for the effect of repeated observations on an individual, random effect models was fitted. \\
\indent To address whether the number of previous episodes affected the incidence of rotaviral diarrhoea, a univariable Poisson regression with number of previous episodes as an ordered categorical variable was performed first.
Then a multivariable model, with previous events as a binary variable was fitted by Poisson regression with random effects.
Records with missing data were dropped from all multivariable analyses.\\
\indent Univariable analyses were performed in \Sexpr{substring(R.version.string,1,16)}, Mantel-Haenszel stratified analyses and multivariable poisson regression were performed in Stata version 12.\\
%
\section{Results}
\indent A total of \Sexpr{length(first$case)} children were recruited to the study, of these \Sexpr{sum(first$case)} (\Sexpr{round((sum(first$case)/length(first$case)*100), 2)} \%) developed rotaviral diarrhoea.
A total of \Sexpr{sum(all$case, na.rm=TRUE)} episodes of rotaviral diarrhoea were recorded over the duration of the study.
The incidence rate for any case of rotaviral diarrhoea was \Sexpr{round((sum(first$case)/(sum(first$fu.time)/365.25))*1000, 2)} per thousand person-years and the incidence rate for all cases was \Sexpr{round((sum(first$max)/sum(first$tyar))*1000 ,2)} per thousand person-years.
\Sexpr{sum(all$case[all$numprev==0])} cases occurred in children who had not previously experienced rotaviral diarrhoea, \Sexpr{sum(all$case[all$numprev==1])} cases occurred in children with one previous episode, \Sexpr{sum(all$case[all$numprev==2])} cases in children with two, \Sexpr{sum(all$case[all$numprev==3])} in children with three, \Sexpr{sum(all$case[all$numprev==4])} in children with four and \Sexpr{sum(all$case[all$numprev==5])} cases in children with five. \\
\indent \Sexpr{length(all$id[all$exittype=="Lost to followup"])} (\Sexpr{round(((length(all$id[all$exittype=="Lost to followup"])/length(unique(all$id)))*100),0)} \%) children were lost to follow-up.
Loss to follow-up was higher in children whose mothers had no education (23 \%) compared to those with primary and middle school education (11 \%) and compared to those whose mothers have higher secondary education (10 \%, p = 0.002).
\Sexpr{length(unique(first$id[is.na(first$lbw)==TRUE]))} children had no data on birthweight.
There was no evidence that the rate of diarrhoea was higher among these children. \\
\indent In univariable analyses there was strong evidence that the rate of diarrhoea was higher in males than in females, was higher in children whose mother had received no education compared to those whose mother had received higher education and was higher for children from households in which beedi were made than those from households in which they weren't (table \ref{epic}). \\
\indent Mantel-Haenzsel analysis indicated interactions between animal ownership and socioeconomic status and between low birthweight and neonatal rotaviral diarrhoea. \\
\indent By univariable analysis there was some evidence that the rate of rotaviral diarrhoea was lower in December 2002 to February 2003 and March to May 2003 compared to September to November 2003 (table \ref{trend}).
After adjusting for subject age, there remained strong evidence that the rate of rotaviral diarrhoea was different from baseline for a number of periods, and that the fit of the multivariable model was better than that of the univariable model. \\
\indent After adjusting for other factors, there was strong evidence that the rate of rotaviral diarrhoea was higher in male children than in female, and higher in households where beedi were made compared to those where beedi were not made (table \ref{multivar}, model 1).
Including interaction terms prevented Stata from estimating maximum likelihood values and so were not included in the final models.
There was strong evidence that having four previous episodes increased the likelihood of experiencing another (table \ref{multivar}, model 2).
There was no evidence that having any other number of previous episodes increased the likelihood of experiencing another compared to having no previous episodes.
Finally, after adjusting for other factors, there was no evidence that the rate of rotaviral diarrhoea was different in children who had previous episodes compared to those who had not (table \ref{multivar}, p = 0.58 lr test of model 1 vs model 3). \\
%
\section{Discussion}
\indent In this study rotaviral diarrhoea was a prevalent disease, affecting over 35 \% of children in the study.
After adjustment for other factors, being of male sex and living in a household where beedi were made were associated with higher rates of disease.
There was strong evidence that the rate of rotaviral diarrhoea varied by time after adjusting for a child's age and there was no evidence that previous episodes of diarrhoea were associated with the probability of a subsequent episode occurring. \\
\indent One limitation is that due to software limitations, it was not possible to study interactions in multivariable models.
This is likely to have concealed some important combinations of conditions, such as being of low birth weight and having neonatal rotaviral diarrhoea, which might have been associated with increased rates of disease. \\
\indent These results may not be reliable if unduly affected by chance, bias or confounding.
Chance is unlikely as there was strong evidence against the null hypothesis for both sex and for beedi manufacture.
Ten variables were tested for association, against a population of 447 children with 236 events. This satisfies the rule of 10 events per regression parameter, suggesting that type I error is unlikely.
The sample size of the study was reasonable and therefore type II error is unlikely. \\
\indent 73 \% of eligible children were recruited into the study.
No reason was given for non-enrollment and thus it is difficult to state whether selection bias occurred.
Had those children who were not enrolled experienced different rates of disease and different exposures, the associations determined here would not be reliable.
Loss to follow-up in this study was higher in children whose mothers had no education, suggesting a possible source of bias.
However, mother's education level was not associated with a different rate of disease and the loss to follow up is not likely to result in bias.
This study employed stool testing to confirm rotaviral diarrhoea.
The method of this testing is not known, but errors in testing could affect the strength of relationships determined.
Random error in diagnosis resulting in under diagnosis would bias any associations towards the null.
Non-random error in diagnosis could result in over or underestimation of associations.\\
\indent Unmeasured confounding could account for some associations observed here.
The reason for the association between male sex and rotaviral diarrhoea is not clear and if another co-mobidity which predisposed to rotaviral diarrhoea were more prevalent among males, then it could account for the difference.
Similarly the mechanism by which beedi manufacture increases the rate is not known but potential major confounders such as SES, education level and household size have been controlled for. \\
\indent In conclusion, the rate of rotaviral diarrhoea varied by time, risk factors included beedi work and male sex but not previous experience of rotaviral diarrhoea.
<<epic_table, results='asis', echo=FALSE>>=
source("epictable.R")
print(tab, include.rownames=FALSE, booktabs=TRUE, size="small", floating.environment='sidewaystable', sanitize.text.function = function(x){x})
rm(tab)
@
\begin{sidewaystable}[p]
\begin{center}
\begin{tabular}{rrrrrrrrr}
\toprule
\multicolumn{1}{l}{} & & \multicolumn{1}{l}{} & & \multicolumn{1}{l}{} & \multicolumn{ 2}{c}{Univariable analysis*} & \multicolumn{ 2}{c}{Multivariable analysis**} \\
\multicolumn{1}{c}{Quarter} & \multicolumn{1}{c}{\shortstack{Median age \\ months (IQR)}} & \multicolumn{1}{c}{\shortstack{Incident \\cases}} & \multicolumn{1}{c}{\shortstack{Person time \\years}} & \shortstack{Rate per \\ 1000 \\person years} & \multicolumn{1}{c}{\shortstack{Rate ratio \\(95 \% CI)}} & \shortstack{Wald \\p-value} & \multicolumn{1}{c}{\shortstack{Rate ratio \\(95 \% CI)}} & \multicolumn{1}{c}{\shortstack{Wald \\p-value}} \\
\midrule
Mar-May 2002 & 1.61 (1.43 - 1.79) & 0 & 0.0011 & 0.00 & 0.00 & 0.999 & 0.00 & 0.999 \\
Jun-Aug 2002 & 2.80 (1.95 - 4.00) & 4 & 0.0131 & 305.78 & 1.67 (0.57 - 4.94) & 0.352 & 1.39 (0.45 - 4.28) & 0.570 \\
Sep-Nov 2002 & 4.32 (2.39 - 6.14) & 9 & 0.0293 & 307.42 & 1.68 (0.76 - 3.74) & 0.203 & 1.33 (0.58 - 3.05) & 0.500 \\
Dec 2002 - Feb 2003 & 4.84 (3.16 - 8.18) & 45 & 0.0536 & 840.17 & 4.60 (2.66 - 7.94) & \textless0.001 & 3.65 (2.06 - 6.49) & \multicolumn{1}{c}{\textless0.001} \\
Mar-May 2003 & 6.87 (4.05 - 10.05) & 38 & 0.0759 & 500.58 & 2.74 (1.56 - 4.80) & \textless0.001 & 2.21 (1.25 - 3.92) & 0.007 \\
Jun-Aug 2003 & 9.07 (5.21 - 12.61) & 39 & 0.0945 & 412.50 & 2.26 (1.29 - 3.95) & 0.004 & 1.97 (1.12 - 3.46) & 0.018 \\
Sep-Nov 2003 & 12.39 (8.68 - 16.00) & 18 & 0.0985 & 182.77 & Reference & & Reference & \multicolumn{1}{l}{ } \\
Dec 2003 - Feb 2004 & 15.55 (11.71 - 19.07) & 26 & 0.0971 & 267.72 & 1.46 (0.80 - 2.67) & 0.213 & 1.90 (1.03 - 3.49) & 0.039 \\
Mar-May 2004 & 18.89 (15.07 - 22.50) & 18 & 0.093 & 193.50 & 1.06 (0.55 - 2.03) & 0.864 & 1.93 (0.97 - 3.83) & 0.060 \\
Jun-Aug 2004 & 21.82 (18.20 - 24.57) & 13 & 0.0809 & 160.67 & 0.88 (0.43 - 1.79) & 0.723 & 2.32 (1.05 - 5.12) & 0.037 \\
Sep-Nov 2004 & 24.18 (21.04 - 26.09) & 12 & 0.0616 & 194.79 & 1.07 (0.51 - 2.21) & 0.864 & 3.91 (1.62 - 9.44) & 0.002 \\
Dec 2004 - Feb 2005 & 25.21 (23.18 - 26.09) & 6 & 0.0379 & 158.46 & 0.87 (0.34 - 2.18) & 0.762 & 3.87 (1.27 - 11.76) & 0.017 \\
Mar-May 2005 & 26.09 (26.09 - 26.09) & 8 & 0.0227 & 351.80 & 1.92 (0.84 - 4.43) & 0.123 & 8.6 (3.04 - 24.31) & \multicolumn{1}{c}{\textless0.001} \\
\bottomrule
\end{tabular}
\end{center}
\caption{Secular trends in rate of rotaviral diarrhoea in children aged two months to two years, India 2002 - 2005. Poisson regression with random effects. IQR: interquartile range. *p for model \textless0.001. **Adjusted for age, p for model \textless 0.001.}
\label{trend}
\end{sidewaystable}
<<trendtable, results='asis', echo=FALSE>>=
#source("incidence_rate.R")
#print(t, include.rownames=FALSE, booktabs=TRUE, sanitize.text.function = function(x){x})
#rm(t)
@
\begin{sidewaystable}[p]
\begin{center}
\begin{tabular}{llllllll}
\toprule
& & \multicolumn{ 2}{c}{Model 1} & \multicolumn{ 2}{c}{Model 2} & \multicolumn{ 2}{c}{Model 3} \\
& & \shortstack{Rate ratio \\(95 \% CI)} & \shortstack{Wald \\p-value} & \shortstack{Rate ratio \\(95 \% CI)} & \shortstack{Wald \\p-value} & \shortstack{Rate ratio \\(95 \% CI)} & \shortstack{Wald \\p-value} \\
\midrule
Sex & Male & Reference & & & & Reference & \\
& Female & 0.68 (0.50 - 0.91) & \multicolumn{1}{r}{0.009} & & & 0.70 (0.52 - 0.93) & \multicolumn{1}{r}{0.015} \\
\shortstack{Number of \\ previous episodes} & \multicolumn{1}{r}{0} & & & Reference & & & \\
& \multicolumn{1}{r}{1} & & & 1.29 (0.95 - 1.75) & \multicolumn{1}{r}{0.104} & & \\
& \multicolumn{1}{r}{2} & & & 0.91 (0.51 - 1.59) & \multicolumn{1}{r}{0.731} & & \\
& \multicolumn{1}{r}{3} & & & 1.47 (0.47 - 4.61) & \multicolumn{1}{r}{0.509} & & \\
& \multicolumn{1}{r}{4} & & & 11.64 (2.89 - 46.97) & \multicolumn{1}{r}{0.001} & & \\
& \multicolumn{1}{r}{5} & & & 0 (0 - 0.93) & \multicolumn{1}{r}{0.985} & & \\
Any previous episode & No & & & & & Reference & \\
& Yes & & & & & 1.26 (0.60 - 2.64) & \multicolumn{1}{r}{0.540} \\
\shortstack{Neonatal rotaviral \\ diarrheoa} & No & Reference & & & & Reference & \\
& Yes & 1.40 (0.70 - 2.79) & \multicolumn{1}{r}{0.344} & & & 1.38 (0.72 - 2.65) & \multicolumn{1}{r}{0.334} \\
Mother's education & None & Reference & & & & Reference & \\
& Primary and middle school & 0.83 (0.57 - 1.21) & \multicolumn{1}{r}{0.340} & & & 0.84 (0.59 - 1.21) & \multicolumn{1}{r}{0.351} \\
& Higher secondary and college & 0.85 (0.59 - 1.23) & \multicolumn{1}{r}{0.392} & & & 0.87 (0.61 - 1.24) & \multicolumn{1}{r}{0.427} \\
Household size & <=5 & Reference & & & & Reference & \\
& >5 & 1.17 (0.86 - 1.59) & \multicolumn{1}{r}{0.326} & & & 1.17 (0.87 - 1.57) & \multicolumn{1}{r}{0.297} \\
Socioeconomic status & Class I & Reference & & & & Reference & \\
& Class II & 0.86 (0.63 - 1.18) & \multicolumn{1}{r}{0.350} & & & 0.87 (0.64 - 1.17) & \multicolumn{1}{r}{0.347} \\
Birthweight & >=2.5 & Reference & & & & Reference & \\
& <2.5 & 0.92 (0.59 - 1.43) & \multicolumn{1}{r}{0.695} & & & 0.92 (0.60 - 1.39) & \multicolumn{1}{r}{0.682} \\
Animal ownership & No & Reference & & & & Reference & \\
& Yes & 0.89 (0.56 - 1.41) & \multicolumn{1}{r}{0.615} & & & 0.90 (0.58 - 1.41) & \multicolumn{1}{r}{0.653} \\
\shortstack{Beedi made \\ in household} & No & Reference & & & & Reference & \\
& Yes & 1.48 (1.08 - 2.01) & \multicolumn{1}{r}{0.014} & & & 1.45 (1.08 - 1.96) & \multicolumn{1}{r}{0.014} \\
\bottomrule
\end{tabular}
\end{center}
\caption{Multivariable poisson regression of potential risk factors for rotaviral diarrhoea in children aged two months to two years, India 2002 to 2005.
Model 1: Poisson regression with random effects, adjusted for age and calendar time, not adjusted for number of previous episodes or any previous episode, p for model < 0.001. p$\theta$ = 0.004.
Model 2: Poisson regression for number of previous episodes, no random effects, adjusted for calendar time and age, p for model 0.071.
Model 3: Poisson regression with random effects, adjusted for age and calendar time, and any previous episode, p for model < 0.001. p$\theta$ = 0.004.p$\theta$ = 0.268. }
\label{multivar}
\end{sidewaystable}
\end{document}