\documentclass[11pt,a4paper]{article}
%\usepackage{palatino}
\usepackage{graphicx}
\usepackage{rotating}
\usepackage{amsmath}
\usepackage{fancyvrb}
\setlength{\parindent}{0pt}
\setlength{\parskip}{6pt}
\setlength{\arrayrulewidth}{1pt}
\renewcommand{\textfraction}{0.05}
\setlength{\textwidth}{16cm}
\setlength{\textheight}{24.5cm}
\setlength{\oddsidemargin}{0pt}
\setlength{\topmargin}{0pt}
\setlength{\headheight}{0pt}
\setlength{\headsep}{0pt}
\newcommand{\bm}[1]{\boldsymbol{#1}}
\newcommand{\var}[1]{\mbox{var}\left({#1}\right)}
\newcommand{\covar}[1]{\mbox{cov}\left({#1}\right)}
\newcommand{\corr}[1]{\mbox{corr}\left({#1}\right)}
\newcommand{\vcite}[2][~]{#1$[#2]$}
%\newcommand{\vcite}[2][]{#1$^{\bm{#2}}$}
\linespread{1.5}
\renewcommand{\baselinestretch}{1.5}

\newcommand{\D}{\displaystyle}  %colocar no numerador e denominador das frações em tabelas para fica com o tamanho normal

\usepackage{fancyvrb}
\DefineVerbatimEnvironment{Rin}{Verbatim}{fontshape=sl}
\DefineVerbatimEnvironment{Rout}{Verbatim}{}
\newcommand{\R}{{\normalfont\textsf{R}}{}}
\newcommand{\code}[1]{\textsl{\texttt{#1}}}
\newcommand{\codeout}[1]{\texttt{#1}}


\title{Comparing modelling strategies for longitudinal data: \\the effect of Vitamin A in childhood diarrhoeal disease}
\author{ V. L. Andreozzi$^\dag$ \& M. S. Carvalho \\[3pt]
{\small \em National School of Public Health - Fiocruz, Brazil.}\\[6pt]
M. L. Barreto\\[3pt]
{\small \em Federal University of Bahia, Brazil.}
}
\date{}

\begin{document}

\maketitle
\renewcommand{\footnoterule}{}
\renewcommand{\thefootnote}{}
\footnotetext{\hspace*{-0.23in}\dag {\em Correspondence to:}
Dr V.L. Andreozzi, National School of Public Health, Fiocruz. Rua Leopoldo Bulh\~{o}es 1480/827, Manguinhos, Rio de Janeiro, RJ, Brazil, Cep: 21041-210. {\tt
E-mail:valeska@fiocruz.br}}


{\linespread{1}
\begin{abstract}

\noindent \textbf{Introduction}:The choice of modelling strategies can improve or impair the best designed studies. It is well known the role of recent developed techniques, such as generalised estimating equations and random effects models, on analysing longitudinal data. On the other hand, late updates on survival analysis models allow the inclusion of multiple time-to-event periods and random effects. In this paper we compare the use of these methods in a longitudinal study on the effects of Vitamin A in preventing childhood diarrhoeal disease. 

\textbf{Methods}: Data comes from a randomised, double-blind, placebo-controlled clinical trial where over 1000 children aged 6-48 months were assigned to receive either vitamin A or placebo, every 4 months over 1 year. Demographic, clinical, nutritional and socio-economic characteristics were collected at the baseline. The subjects were visited at home three times a week, to register the number of liquid or semi-liquid motions a day. Three response variables were considered: number of liquid or semi-liquid motions per week, numbers of days with diarrhoea per week, in the longitudinal models; numbers of days without diarrhoea in the survival models. 

% \textbf{Results}: The survival models are more comparable to the longitudinal models for the number of days with diarrhoea, since both of them model the occurrence of diarrhoea episodes. For the two ways of looking at the problem, no significant effect was found for the vitamin A supplementation. The only significant effect was found when modelling the severity of the disease, using the number of liquid or semi-liquid motions per week. The most important covariates were the age of the child, as expected, and the number of domestic appliances (refrigerator, radio, TV). 

\textbf{Conclusion}: The comparison of the models showed reasonable consistency; even the differences found can be attributed to different measures and ways of looking at the problem. The major benefit trying those various models is the ability of detecting small change in point and confidence interval estimates, including independent variable of interest and all potential risk factors. The variety of statistical models is increasing sharply with the development of statistical tools and availability of software. Unfortunately the choice among them is usually more dependent on previous knowledge then on the real question addressed by the study. Quite often data like this one is only analysed using simple summary measures, in classical models easily adjusted in common statistical packages, causing loss of information so hardly gathered. There is a strong need in the epidemiology environment of incorporation of an entire set of statistic techniques dealing with dependent data, particularly from longitudinal studies. 

\paragraph{Keywords: Longitudinal Study, Marginal models, Random effects models, Survival Models}
\end{abstract}
}



\section{Results} \label{sec:results}

Testing the inclusion of many available covariates, the total time of breastfeeding and several variables on water quality and sewerage systems were not significant in any model. This was due to their small variability: environment conditions were very similar for all households, and the total time of breastfeeding lasted less then 3 months for more then 85\% of the chidren. Among the socio-economic covariates, the number of domestic appliances (variable name: appliance) was associated with the outcome in all models hence being chosen. 

The results of the repeated measures models are presented in Table~\ref{tab:longit}. The point estimates of the marginal models indicates that the odds of a day with diarrhoea increased for the children not taking the vitamin supplement between 14\% (not significant for 5\% CI) and 20\%. The number of domestic appliances was significant in all models, and with equal point estimate, if considering a precision of two digits: the odds ratio for each appliance added to the housedhold decreased the day probability of diarrhoea in 15\%. Age was protective as well, decreasing the odds between 5\% and 6\% per month. Comparing the marginal models, we observe that, in general, the auto-regressive correlation structure model is slightly more conservative then the exchangeable one, with both the Vitamine A effect and linear calendar time not significant within 95\% of confidence. The loss of significance of time in the AR model is due to this structure, as it, in part, model the slow decrease of the probability of days with diarrhoea. Auto-regressive models apprehend, in part, the time trend, which is the reason for differencing time series before fitting AR models. And that is the reason as well to give some preference to this model over the EX one, as it controls for the time variation. 

Both marginal models were more conservative then the random effects one. Adjusting he random effects model with non-linear time trend included the mother's shooling, increasing the odds as the schooling years decreased. The density of the posterior distributions, a very informative result allowed by Bayesian inference, are presented in Figure \ref{fig:densit}. The time trend, presented in Figure \ref{fig:trend}, is non-linear, and similar to the explo one??????? 
% cade o grafico?  

Using the survival approach (Table~\label{tab:surv}), the measure is now time-to-event, a diarrhoa episode. The measure of association estimated in the model ($exp(\beta)$) is the hazard ratio. Vitamin A was marginally significant, with a very small effect: approximately 8\% -- 9\% of diarrhoea episodes prevented. The effect of age was smaller as well, decreasing the hazard of a new episode of diarrhoea in 3\% per month. Again, the number of domestic appliances was a stronger predictor, and the mother education was significant only for the smallest number of formal education as compared with high attendance. The main difference between the marginal and the frailty models is the proportion of variability explained by the model, estimated in just 10\% for the marginal model and 53\% for the frailty one. Time could not be included in the model, because the use of a counting process model, important to deal with time varying covariates such as age defines each event in a time point instead of using the more common total time of observation. Thus, the happening is exactly at the same moment then the calendar time covariate would be, hence unidentifiable. 





\begin{sidewaystable} \label{tab:longit}
\begin{center}
\caption{\small Marginal models with different correlation structure and random effects models estimates. (OR=Odds Ratio, SE=Standard Error, SD=Standard Deviation,$\phi$ = dispersion parameter,$\alpha$=correlation parameter, $\sigma_{\gamma_i}^2$ =intercept variance )}
\begin{tabular}{l c c |  c c  |  c c}
\hline\hline  
 & \multicolumn{2}{c|}{Marginal - Exchangeable} & \multicolumn{2}{c|}{Marginal - Auto-regressive} & \multicolumn{2}{c}{Random Effects} \\ \hline
% & OR &  95\% Confidence Interval &   Estimate & 95\% Confidence Interval &   OR & 95\% Credible Interval \\ \hline 
  & OR &  95\% Confidence   &    OR & 95\% Confidence  & OR & 95\% Credible  \\ 
             &  &   Interval &      &  Interval &  &  Interval \\ 
            \hline 
Intercept    & 0.307  &  [0.232; 0.407]    &      0.285  & [0.214; 0.380]  &  0.185  &  [0.136 ; 0.252]   \\
Group (Vit)  & 0.828  &  [0.708; 0.968]    &      0.876  & [0.749; 1.024]  &  0.847  &  [0.720 ; 0.997]   \\
Sex (M)      & 1.061  &  [0.909; 1.239]    &      1.102  & [0.945; 1.284]  &  1.057  &  [0.898 ; 1.244]   \\
Age (months) & 0.961  &  [0.955; 0.968]    &      0.959  & [0.953; 0.966]  &  0.951  &  [0.945 ; 0.957]   \\
Time (week)  & 0.996  &  [0.992; 0.999]    &      0.998  & [0.995; 1.002]  &  \multicolumn{2}{c}{non-linear} \\
Schoolmom2   & 1.045  &  [0.860; 1.269]    &      1.086  & [0.898; 1.314]  &  1.241  &  [1.017 ; 1.515]   \\
Schoolmom3   & 1.152  &  [0.904; 1.467]    &      1.247  & [0.968; 1.606]  &  1.327  &  [1.040 ; 1.695]   \\
Appliance    & 0.849  &  [0.795; 0.907]    &      0.852  & [0.797; 0.911]  &  0.845  &  [0.784 ; 0.911]   \\ \hline 
& \multicolumn{2}{l|}{$\phi$= 0.536 (SE=0.080)} & \multicolumn{2}{l|}{$\phi$= 0.578 (SE=0.072)} & 
	\multicolumn{2}{l}{$\sigma_{\gamma_i}^2$ =1.663 (SD=0.087)}\\
& \multicolumn{2}{l|}{$\alpha$= 0.146 (SE=0.032)} & \multicolumn{2}{l|}{$\alpha$= 0.527(SE=0.055)} \\ 

%& \multicolumn{2}{l|}{Dispersion ($\phi$): 0.536 (SE=0.080)} & \multicolumn{2}{l|}{Dispersion ($\phi$): 0.578 (SE=0.072)} & 
%	\multicolumn{2}{l}{Intercept variability($\sigma_{\gamma_i}^2$) =}\\
%& \multicolumn{2}{l|}{Correlation ($\alpha$): 0.146 (SE=0.032)} & \multicolumn{2}{l|}{Correlation ($\alpha$): 0.527(SE=0.055)} \\ 

\hline
\end{tabular}      
\end{center}    
\end{sidewaystable}


\begin{sidewaystable} \label{tab:surv}
\begin{center}
\caption{\small AG and Frailty survival models estimates (HR= hazard ratio, CI= Confidence Interval, p=p-value) }
\begin{tabular}{l c c  | c c }
\hline\hline  
 & \multicolumn{2}{c|}{AG Model}            & \multicolumn{2}{c}{Frailty Model} \\ \hline                        
              & HR &  95\% CI   &               HR &  95\% CI     \\ \hline 
Group (Vit)   &  0.915 &  [0.817; 1.025]       &  0.914  &  [0.809; 1.03]      \\         
Sex (M)       &  1.039 &  [0.928; 1.164]       &  1.033  &  [0.914; 1.17]      \\         
Age           &  0.973 &  [0.968; 0.978]       &  0.975  &  [0.970; 0.98]      \\         
Schoolmom2    &  1.121 &  [0.969; 1.296]       &  1.113  &  [0.958; 1.29]      \\         
Schoolmom3    &  1.206 &  [1.010; 1.439]       &  1.222  &  [1.019; 1.47]      \\         
Appliance     &  0.896 &  [0.851; 0.942]       &  0.886  &  [0.836; 0.94]      \\ \hline
       
% &  & & \multicolumn{2}{l}{ p-value - Frailty(child)  0.000}    \\ 
 & \multicolumn{2}{l|}{R$^2$ = 0.104}              & \multicolumn{2}{l}{R$^2$ = 0.528 (Maximum=1)}\\
 &   &  & \multicolumn{2}{|l}{Frailty Variance= 0.836 (p$<$0.001)}\\
% & \multicolumn{2}{l|}{LR test= 865 (df=6, p=0)}     &  \multicolumn{2}{l}{LR test= 5939 (df=881 , p=0)}\\
% & \multicolumn{2}{l|}{Wald test = 178 (df=6 , p=0)} &  \multicolumn{2}{l}{Wald test = 148 (df=881,p=1)}\\          
% & \multicolumn{2}{l|}{Score test = 870 (df=6, p=0), Robust = 129 p=0}   \\ 
\hline    
\end{tabular}      
\end{center}
\end{sidewaystable}






\begin{figure} \label{fig:densit}
\begin{center}
\includegraphics[angle=-90,width=6cm]{m15densgroup.eps} 
\includegraphics[angle=-90,width=6cm]{m15denssex.eps}
\includegraphics[angle=-90,width=6cm]{m15denssch2.eps}
\includegraphics[angle=-90,width=6cm]{m15denssch3.eps}
\includegraphics[angle=-90,width=6cm]{m15densage.eps}
\includegraphics[angle=-90,width=6cm]{m15densapli.eps}
\caption{\small Probability Density function of REM coefficients. Center line is the estimated $\bm{\beta}$ and the dotted segment is the 95\% credible interval. Black dot on horizontal axis calls attention to the $\beta$ = 0.}
\end{center}
\end{figure}


\begin{figure} \label{{fig:trend}}
\begin{center}
\includegraphics[angle=-90,width=10cm]{m15week.eps}
\caption{\small Non-linear Effect of time. Shown are the posterior means together with 95\% and 80\% pointwise credible intervals.}
\end{center}
\end{figure}







\section{Discussion} \label{sec:discussion}

Recent published paper (ref Tooth) presented a thorough systematic review on the quality of reporting of observational longitudinal research, where although 98\% of the selected papers mentioned the specific statistical methods by name, almost 25\% did not stated the use of specific methods to analyse longitudinal studies. And this aspect was not debated as carefully as other aspects. Our paper deals with this specific point: not only 25\% of reports ignores specific aspects of appropriate "longitudinal" methods, but such an useful article ignores almost ignores this problem as well. It was our intent to bring an applied example of the variety of potencially useful methods, contributing to bridge the gap between statisticians and epidemiolgists.

Instead of dealing with longitudinal information as just as a "treated" versus "non treated", or "present" versus "absent" characteristic, the conceptualization and statistical modeling of change in the form of a trajectory as a function of time offer more efficient strategies to understand the effect of any covariate, either in intervention or observational longitudinal studies. 




The first part of this formulation is similar to that of the standard GLM, but the important distinction is that it involves the \textit{conditional} mean and variance of the response given values for the random effects. This may be contrasted with the GLMs in Section~\ref{subsec:marginal} which model the \textit{marginal} mean and variance of the response i.e. that  averaged over the probability distribution for any unobserved random effects (subject heterogeneity). In the random effects model a within-subject covariance structure is indirectly generated by assuming that each individual has unique random effects, distributed around the fixed population effects\vcite{20}. In the marginal model subject-specific random effects are `averaged out' of the model, within-subject correlation being handled via a function specified separately to the mean structure. The important point here is that there is a subtle difference in the interpretation of fixed effect parameters in the marginal and random effect models. In the former, they are `population averaged' effects applicable to groups of subjects exhibiting similar heterogeneity to that observed amongst the subjects in the study. In the random effects model the fixed parameters are interpreted as effects on an
`average subject' in the study. Some have argued that the latter inference is more useful in clinical trials, where assessment of treatment effect on a average patient is of more interest than its average effect on a population as a whole\vcite{21}.


All this modelling exercise allows us to make a few assertions:
\begin{itemize}
\item Age is by far the most important covariate in diarrhoeal disease determination.
\item The best family socio-economic indicator is the number of household appliances.
\item The mother's education level is not as associated with diarrhoeal disease in this community as described in the literature, in spite of a significant number of mothers in all categories (281 with no education, 588 with primary school and 326 with high school). 
\item Calendar time, or time since the beginning of the study, when modelled, presented a strong protective effect. It should be noticed that it was included in the models in a linear way, and visual inspection of the trend showed a non linearity.
\end{itemize}

A significant effect of vitamin A could be only detected preventing the total number of liquid and semi-liquid motions, and not for days with the disease neither for days free from the disease. Considering the severity models, the random effects one presents one advantage over marginal model: it is able to capture the larger variability among children, via the random intercept. When marginal model is used the interpretation of the coefficients combines the within-subject relationship with a between-subjects relationship.

The model for days with diarrhoea presented different effects for our main interest, depending on the working correlation structure: significant effect when exchangeable structure were adopt, and not significant when an auto-regressive was chosen. Unfortunately there is not a specific method to choose the correlation structure. When the responsible variable is continuous, the correlation structure of the observed data serves as a guide for its choice, not our case. Considering the literature our results are moving in the opposite direction because the GEE analysis was not robust against the wrong choice of a correlation structure for non-continuous data. So the choice of the correlation must be made based on the knowledge of the researcher on problem to be analized.

The most interesting feature of the survival model with random effects is its ability to capture and measure through the R$^2$, the unexplained variability, in spite of its limit in including the calendar time. It should be noted that 
usually at this point the question that arises is: which is the best model? How to choose among them? First of all we should discuss in theory which one is a better approach to the problem. Besides, many improvements could be made to those models using the same data set, for instance, the inclusion of the time in a non linear way, as well as to the study design, in order to treat the possible learning effect detected via the decay time trend. What we want to emphasise is that \textit{All models are wrong. Some models are useful} \cite{box79}. And what we really want to know is: should we propose the supplementation of vitamin A as an strategy for controlling diarrhoeal disease in deprived populations? Well, the answer is not as straightforward as we would like. Although we have some evidence in favour of a real protective effect of vitamin A on the severity of diarrhoea, even when controlled by socio-economic level, mother's education and so on, the time trend can be a strong confounder, and should be properly addressed. Besides, obviously a cost-effectiveness study should be done, comparing the efficacy of this strategy with alternatives such education, casually detected via the time trend. Statistical modelling is a tool in medical sciences and an art. No model itself is right, and we should consider previous knowledge, robustness, performance using different suppositions, instead of just looking at p-values and positive results.

\bibliography{tese}
\bibliographystyle{unsrt}


\end{document}
