\documentclass[a4paper, 10pt]{report} \usepackage{natbib} \bibliographystyle{plainnat} \usepackage[T1]{fontenc} \usepackage{url} \usepackage[utf8]{inputenc} \usepackage{graphicx} %\usepackage{tikz} %\usetikzlibrary{decorations,arrows,shapes} \usepackage[margin=0.9in]{geometry} %\usepackage[left=3cm,right=3cm,top=2cm,bottom=2cm]{geometry} \usepackage{hyperref} \usepackage{listings} \usepackage{xspace} \usepackage{tabularx} \usepackage{makeidx}\makeindex \usepackage{algorithmic} \usepackage{algorithm} \usepackage{amsmath,amsthm,amsfonts,amssymb} \setlength{\parindent}{0mm} \setlength{\parskip}{1mm} \newcommand{\commentout}[1]{} \renewcommand{\theequation}{\thesection.\arabic{\equation}} \numberwithin{equation}{section} \theoremstyle{definition} \newtheorem{Def}{Definition}[section] \newtheorem{Rem}[Def]{Remark} \newtheorem{RemDef}[Def]{Remark und Definition} \newtheorem{DefRem}[Def]{Definition und Remark} \newtheorem{Example}[Def]{Example} \theoremstyle{plain} \newtheorem{Theorem}[Def]{Theorem} \newtheorem{DefTheorem}[Def]{Definition and Theorem} \newtheorem{Corollary}[Def]{Corollary} \newtheorem{Lemma}[Def]{Lemma} \newcommand{\R}{\ensuremath{\mathbb{R}}\xspace} \newcommand{\NN}{\ensuremath{\mathbb{N}_0}\xspace} \newcommand{\N}{\ensuremath{\mathbb{N}}\xspace} \newcommand{\sF}{\ensuremath{\mathcal{F}}\xspace} \newcommand{\sN}{\ensuremath{\mathcal{N}}\xspace} \newcommand{\Pot}{\ensuremath{\mathfrak{Pot}}\xspace} \newcommand{\kronecker}{\raisebox{1pt}{\ensuremath{\:\otimes\:}}} \DeclareMathOperator{\range}{range} \DeclareMathOperator{\Var}{Var} \DeclareMathOperator{\trace}{trace} \newcommand{\skp}[1]{\left\langle#1\right\rangle} \renewcommand{\epsilon}{\varepsilon} \renewcommand{\phi}{\varphi} \newcommand{\id}{\text{id}} \newenvironment{Proof}{\par\noindent\upshape\textit{Proof. }\nopagebreak}{\qed\par} \usepackage{setspace} \onehalfspacing \begin{document} %\VignetteEngine{knitr::knitr} %\VignetteIndexEntry{Crossover - A search algorithm and GUI for cross-over designs} \title{Crossover - A search algorithm and GUI for cross-over designs} \author{Kornelius Rohmeyer} \maketitle %\newpage \tableofcontents <>= # knitr has to be loaded for 'set_parent' and CRAN checks and also for opt_chunk during build process. library(knitr) if (exists("opts_chunk")) { opts_chunk$set(concordance=TRUE) opts_chunk$set(tidy.opts=list(keep.blank.line=FALSE, width.cutoff=80)) #opts_chunk$set(size="footnotesize") #opts_chunk$set(size="tiny") opts_chunk$set(size="scriptsize") opts_chunk$set(cache=TRUE) opts_chunk$set(autodep=TRUE) } # See http://yihui.name/knitr/hooks for the following code. CrossoverNamespace <- function(before, options, envir) { if (before) { ## code to be run before a chunk #attach(loadNamespace("Crossover"), name="namespace:Crossover", pos=3) attach(list2env(as.list(asNamespace("Crossover"))), name="exposed:Crossover", pos=3, warn.conflicts=FALSE) } else { ## code to be run after a chunk detach("exposed:Crossover") } } if (exists("opts_chunk")) { knit_hooks$set(withNameSpace = CrossoverNamespace) } library(Crossover, quietly=TRUE) options(width=80) options(digits=4) startGUI <- function(...) {invisible(NULL)} #options(prompt="> ", continue="+ ") library(MASS) library(multcomp) library(ggplot2) library(Matrix) library(Rcpp) bibCall <- TRUE @ \newpage \chapter{Introduction} This package provides more than two hundred cross-over design from literature, a search algorithm to find efficient cross-over designs for various models and a graphical user interface (GUI) to find/generate appropriate designs. The computationally intensive parts of the package, i.e.\ the search algorithm, is written using the R packages Rcpp and RcppArmadillo (\cite{eddelbuettel2011rcpp} and \cite{rcpparmadillo}). The GUI is written in Java and uses package rJava (\cite{rJava}). \section{Installation} %If you don't already have R on your system, you can download a bundled %version of R and Crossover from \url{http://www.algorithm-forge.com/gMCP/bundle/}. %Open R and type \texttt{install.packages("Crossover")}, %select an arbitrary mirror and Crossover will be downloaded and installed. Once it is installed, whenever you start R you can load the Crossover package by entering \texttt{library(Crossover)} into the R Console. The graphical user interface as shown in figure \ref{GUI} is started with the command \texttt{CrossoverGUI()}. \begin{figure}[ht] \centering \includegraphics[width=0.8\textwidth]{figures/GUI.png} \caption{\label{GUI} Cross-Over Design GUI.} \end{figure} %If you run into problems, see %\url{http://cran.r-project.org/web/packages/gMCP/INSTALL} or please write us %an email at %\href{mailto:help@small-projects.de}{\texttt{help@small-projects.de}}. %We are eager to help and to learn about existing problems. \section{Overview} The catalogue, collected and compiled by Professor Byron Jones, contains 241 designs from the following literature: \cite{anderson2002locally}, \cite{archdeacon1980some}, \cite{atkinson1966designs}, \cite{balaam1968two}, \cite{berenblut1964change}, \cite{blaisdell1980partially}, \cite{davis1969cyclic}, \cite{federer1964tied}, \cite{fletcher1987new}, \cite{iqbal1994efficient}, \cite{lewis1988factorial}, \cite{cochran1941double}, \cite{patterson1962change}, \cite{pigeon1985residual}, \cite{prescott1999construction}, \cite{quenouille1953design}, \cite{russell1991construction}, \cite{lucas1956switchback}, \cite{williams1949experimental}, \cite{prescott1994}, \cite{bate2002} Further \Sexpr{length(load(paste(system.file("data", package="Crossover"), "pbib2combine.rda", sep="/")))} designs are constructed from partially balanced incomplete block designs from \cite{clatworthy1973tables} and balanced crossover designs. \subsection{Design selection in the GUI} The GUI will show appropriate designs from the catalogue according to the number of treatments, periods and the range of sequences the user enters.\footnote{A table referencing all available designs and the respective number of treatments, periods and sequences is available by calling \texttt{buildSummaryTable()}.} Further functions from package crossdes (\cite{crossdes}) are called to create designs for the specified values if possible. In figure \ref{GUI} you can see the following four checkboxes, that allow you to see only specific subsets: \begin{description} \item[Designs from package archive] The previously noted designs from literature are shown. \item[Designs generated by package crossdes] Activating this option will result in short delays when displaying the catalogue, since the crossdes algorithms are called. \item[Designs manually entered] All designs entered on tab \emph{"Input own design"} are shown. \item[Designs from previous search runs] All designs from previous search runs are shown. \end{description} The column \text{av.eff.trt.pair} shows the \emph{average effiency of pairwise treatment parameter differences}. To define these efficiencies we compare a design $D$ with an \emph{"ideal design"} $D_I$ without period, subject or carry-over effects, where each treatment $i\in\{1,\ldots,v\}$ occurs $n_i$ times in $D$ and $D_I$. For $D_I$ the variance of the treatment parameter difference estimate $\hat\tau_i-\hat\tau_j$ for treatment $i$ and $j$ is $1/n_i+1/n_j$. We define the efficiency for a design $D$ as the quotient of the variances under $D_I$ and $D$, i.e. \[\frac{\Var_{D_I}(\hat\tau_i-\hat\tau_j)}{\Var_D(\hat\tau_i-\hat\tau_j)} =\frac{\frac1n_i+\frac1n_j}{\Var_D(\hat\tau_i-\hat\tau_j)}.\] %For a design $D$ the efficiency of the estimate of a treatment parameter difference is the quotient of the variances of this The tab \emph{"Input own design"} provides you with the possibility to analyse and save your own designs easily or to use them as starting points for the search algorithm. The drop-down menu for \emph{model} let you specify which model you are interested. These models are described in detail in chapter \ref{chap:models}. In case of the placebo or proportionality model you can specify \emph{further model parameters} (namely the number of placebos and the proportionality parameter, respectively). \newpage \subsection{Designs in R} Our crossover designs in R are numeric matrices, where the elements represent the treatments, the rows represent the periods and the columns the subjects. A data frame referencing all available designs and the respective number of treatments, periods and sequences is available by calling the \texttt{buildSummaryTable()} function. Designs referenced in this table can be accessed via the \texttt{getDesign()} function. For example the Williams design for three treatments is represented by a $3\times6$-matrix and assigns each of the three treatments once to each of the six subjects: <>= getDesign("williams3t") @ Each treatment occurs six times, two times in the each period and each treatment follows each other treatment exactly two times. %\[Y_{ijk}=\mu+\pi_j+\tau_{d[i,j]}+\lambda_{d[i,j-1]}+s_{ik}+e_{ijk}\] If we are interested in the variance of the treatment parameter difference estimates, we can use the function \texttt{general.carryover}: <>= design <- getDesign("williams3t") general.carryover(design, model=1) @ We see that the Williams design is a balanced design. The following nine models which are discussed in chapter \ref{chap:models} are implemented: <>= cat(paste(1:9, ": \"", Crossover:::models, "\"", sep=""), sep="\n") @ \iffalse \subsubsection{Efficiency remarks when selecting a design} Efficiencies for a design $D$ are shown in comparision to an ideal design $D_I$ where each treatment occurs equally often in $D$ and $D_I$. If for example we take a look at the design \texttt{pidgeon1}, we see that treatment 1 occurs nine times, while all other treatments occur only six times: <>= getDesign("pidgeon1") @ The variances of the estimates for the treatment differences are given by: <>= design <- getDesign("pidgeon1") design.efficiency(design)$var.trt.pair.adj @ In an ideal design the variances are: <>= v <- length(levels(as.factor(design))) im <- matrix(0, v, v) vn <- sapply(1:v, function(x) {sum(design==x)}) for (i in 1:v) { for (j in 1:v) { if (i!=j) { im[i, j] <- 1/vn[i]+1/vn[j] } } } im @ We see that the variances for comparing the first treatment to any other treatment are lower (0.3720) than the variances of treatment difference estimates within the last three treatments. <>= design.efficiency(design)$eff.trt.pair.adj @ The average efficiency is given by \Sexpr{mean(design.efficiency(design)$var.trt.pair.adj[upper.tri(design.efficiency(design)$var.trt.pair.adj)])}, \begin{enumerate} \item A higher average efficiency does obviously not mean that a design is better suited for a sepific problem. If you are interested in Dunnett type contrasts you should use other designs as if you are interested in Tukey or Williams type contrasts. \item When comparing two designs with different \item It is better to look directly at the variances of the \end{enumerate} \fi \newpage \subsection{Algorithm Search} \begin{figure}[ht] \centering \includegraphics[width=0.5\textwidth]{figures/Algorithm.png} \caption{\label{Algorithm} Panel for algorthmic search of cross-over designs.} \end{figure} In figure \ref{Algorithm} the preliminary graphical interface for the search algorithm is shown with the following options: \begin{itemize} \item The \emph{covariance pattern} can be \emph{"Independence"}, \emph{"Autoregressive Error"} or \emph{"Autoregressive Error"}. Except for \emph{"Independence"} a \emph{covariance pattern coefficient} $\rho$ has to be specified. %TODO %TODO \item If you \emph{include fixed subject effects in the design matrix} $\ldots$ \item Specify the exact \emph{number of sequences}. (The number of treatments and periods is already specified in the top panel of the GUI.) \item Optionally specify the \emph{exact number of treatment assignments}. The GUI default is to let the algorithm figure out good/optimal assignments. But depending on further information (information from theoretical results or treatments more important than the others\footnote{Different weights of treatment importance should be specified as weighted contrasts. See item \emph{contrasts}.}, etc.) the number of treatment assignments can be specified. \item You can specify that the design should be constructed in a way that \emph{in each sequence/period a treatment occur as evenly as possible}. This restriction will normally decrease the efficiency of the algorithm. \item The GUI default is an all-pair comparison of all treatments with equal weights. Change the \emph{contrast weights} accordingly if you are interested in other contrasts or different weights. %\emph{User defined contrasts} can be used and a R matrix of $n$ %contrasts (i.e. of dimension $n\times v$) will be used. \item Pressing the \emph{"Compute Design"} button will start the search algorithm described in section \ref{sec:search}. After a few seconds the result will be shown in the previous empty text area on the right. \end{itemize} <>= @ <>= @ \iffalse \section{Creation of Crossover Designs from Partially Balanced Incomplete Block Designs with Two Associate Classes (PBIB(2))} We start with a non-crossover PBIB(2) design. A comprehensive source for these is \cite{clatworthy1973tables}. \fi \section{Random Subject Effects Model} See \cite{jones2003design}, 5.3, page 213ff. The model stays the same \[Y_{ijk}=\mu+\pi_j+\tau_{d[i,j]}+\lambda_{d[i,j-1]}+s_{ik}+e_{ijk},\] but we also assume that the subject effects follow a normal distribution: $s_{ik}\sim\sN(0,\sigma_s^2)$. In matrix notation we have \[Y=X\beta+Z\gamma+\epsilon\] with $X$ and $Z$ the fixed and random effects design matrices\footnote{Note that $X$ and $Z$ are different from the ...}, $\epsilon\sim\sN(0,\Sigma)$ and $\gamma\sim\sN(0,D)$. Then\footnote{See for example \cite[chapter 5]{lee2006generalized}.} \[\Var(Y)=ZDZ^T+\Sigma.\] For known $V:=\Var(Y)=ZDZ^T+\Sigma$ the MLE and BLUE is given by \[\hat\beta=(X^tV^{-1}X)^{-1}X^tV^{-1}Y.\] <>= @ \nocite{*} \newpage \addcontentsline{toc}{section}{Index} \printindex %\newpage %\listofalgorithms %\listoffigures %\listoftables \newpage \addcontentsline{toc}{section}{Literatur} \bibliography{literatur} \addcontentsline{toc}{section}{Table of Symbols} \section*{Table of Symbols}\footnotesize %\twocolumn[\section{Symbolverzeichnis}] \begin{tabularx}{\textwidth}{lX} \multicolumn{2}{l}{\textbf{Sets}}\\ \R& set of real numbers\\ \NN& set of natural numbers (including 0)\\ $\Pot(X)$ & power set of set $X$, i.e.\ the set of all subsets of $X$\\ \\ \multicolumn{2}{l}{\textbf{Variables}}\\ $v$& number of treatments\\ $p$& number of periods\\ $s$& number of sequences\\ $\mu$& intercept\\ $\pi_j$& period effect for period $j$\\ $\tau_{d[i,j]}$& direct treatment effect for treatment $d[i,j]$ in period $j$ of sequence $i$\\ $\lambda_{d[i,j-1]}$& first-order carry-over effect (0 for $j-1$=0)\\ $s_{ik}$& $k$th subject effect on sequence $i$\\ $e_{ijk}$& random error with zero mean and variance $\sigma^2$\\ \\ %\end{tabularx}\\ %\begin{tabularx}{\textwidth}{lX} \multicolumn{2}{l}{\textbf{Functions}}\\ $X'$ & transpose of matrix $X$\\ $X^+$ & Moore-Penrose pseudoinverse of $X$\\ $\skp{\cdot,\cdot}$ & standard direct product $\skp{x,y}=\sum_{j=1}^nx_j\cdot y_j$ for $x,y\in\R^n$\\ $\id_X$ & identity on $X$, i.e.\ $\id_X:\;X\rightarrow X,\;x\mapsto x$\\ \\ \multicolumn{2}{l}{\textbf{Other Symbols}}\\ $\sN(\mu,\sigma^2)$ & Normal distribution with mean $\mu$ and variance $\sigma^2$.\\ $\sN(\mu,\Sigma)$ & Multivariate normal distribution with mean $\mu$ and covariance matrix $\Sigma$.\\ \end{tabularx} %\onecolumn \normalsize \newpage \end{document}