%\VignetteIndexEntry{Introduction to R Lab}
%\VignettePackage{IWB2011}

\documentclass{article}

\usepackage{BioconductorVignette}
\SweaveOpts{keep.source=TRUE,width=9,height=8,eps=FALSE,include=FALSE}

\begin{document}

\title{Introduction to R}
\author{Nishant Gopalakrishnan, Martin Morgan}
\date{19-21 January, 2011}
\maketitle
\tableofcontents

<<setup, echo=FALSE>>=
options(width=60)
olocale = Sys.setlocale(locale="C")
vers <- R.Version()
vers <- paste(vers$major, vers$minor, sep = ".")
@

\section{Introduction}

This lab introduces basic \R{} operations by inputing and manipulating
data describing a microarray experiment invovling 128 individuals with
acute lymphoblastic leukemia. Covariates include measures such as age,
sex, type, stage of the disease, etc., and are provided as a comma
separated file \Robject{pData.csv}.  Our goal in this exercise is to
\begin{itemize}
\item Read the covariates into \R{}.
\item Perform data manipulations such as tabulation and subsetting 
\item Visualize some of the data using the \Rpackage{lattice}
  package.
\end{itemize}

\section{Loading tabular data into \R{}}

Here we load the microarray experiment covariates from a `csv'
(comma-separated value) file. The file is located in the
\Robject{extdata} folder of the \Rpackage{IWB2011} package.

\R{} provides several functions such as \Rfunction{read.table} for
reading in meta data files into a \Robject{data.frame} with
appropriate column and row names from the header information provided
in the file. Convenience functions such as \Rfunction{count.fields}
are also available for discovering problems in files (such as certain
rows in the file having different number of fields) when using the
\Rfunction{read.table} function.

\begin{Exercise}
  \begin{itemize}
  \item Start \R{}, and load the \Rpackage{IWB2011} pacakge.
  \item Use the \Rfunction{system.file} to locate the path to the
    files \Robject{exprsMat.csv} and \Robject{pData.csv}
  \item Make use of the \Rfunction{count.fields} function on the
    \Robject{pData.csv} file to ensure that the file has the same
    number of fields in each line of the file.
  \item Use the \Rfunction{read.table} function to read in the
    experimental meta data into an \R{} variable \Robject{pdOrig}.
  \item View the first few records of the data using \Rfunction{head}.
  \item Obtain a brief summary of the data using \Rfunction{summary}.
  \item Tabulate the number of males and females in the study by
    selecting the \Robject{sex} column and using \Rfunction{table}.
  \end{itemize}
\end{Exercise}

\begin{Solution}
<<readData>>=
library(IWB2011)
phenoPath <- system.file( "extdata", "pData.csv", package="IWB2011") 
pdOrig <- read.table(phenoPath)
names(pdOrig)
head(pdOrig)
summary(pdOrig)
table(pdOrig$sex)
@
\end{Solution}

\section{Subset data}

The \Robject{pdOrig} variable is a \Robject{data.frame} with columns
representing the various covariates that describe the experiment (age,
sex, etc.). The row names correspond to the sample Id's. The column BT
is a factor indicating the tumour cell type (B1, B2, T1 ,T2 etc. with
B indicating B-cell type and T indicating T-cell type).  Similarly the
column \Robject{mol biol} is a factor indicating the molecular biology
of the cancer. (BCR/ABL, NEG, E2A/PBX1 etc.) The \Robject{mol biol}
column can be accessed using \Rfunction{pdOrig[["mol biol"]]}.

In this section, we make use of indexing, subsetting, factors etc.  to
modify \Robject{pdOrig} to select a subset of samples.  We are
specifically interested in B-cell tumours with molecular biology type
"NEG" or "BCR/ABL".

\begin{Exercise}
  \begin{itemize}
  \item Identify the samples that are "NEG" or "BCR/ABL" molecular
    biology type using the column \Robject{mol biol}, the
    \Rfunction{\%in\%} function and the \Rfunction{which} functions in
    \R{}.

  \item Use the \Rfunction{grep} function on the \Robject{BT} column in
    the \Robject{pdOrig} \Robject{data.frame} to identify the B cell
    tumours.

  \item Identify samples that are both B cell tumours and are
    "BCR/ABL" or "NEG" using the \Rfunction{intersect} function on the
    indices that we have previously computed.

  \item Subset the phenotypic data \Robject{pdOrig} to create a new
    variable \Robject{psubData}. Note:The rows of the data.frame
    represent samples.
  \end{itemize}
\end{Exercise}

\begin{Solution}
<<subsetPhenoData>>=
types <- c("NEG", "BCR/ABL")
moltyp <- which(as.character(pdOrig$mol.biol) %in% types)
bcell <- grep("^B", as.character(pdOrig$BT))
indx <- intersect(bcell, moltyp)
psubData <- pdOrig[indx,]
@
\end{Solution}
%% $

\section{Recodig factor levels}

The covariate data for some variables, for example \Robject{BT}, is
represented using a variable of type factor. The `levels' are the
distinct categorical values of a factor. Subsetting a factor leaves
the levels of the variable unchanged. In this exercise, we take a look
at the levels of the factor variables that we have just subsetted, and
then update the levels.


\begin{Exercise}
  \begin{itemize}
  \item Observe the levels for the \Robject{mol.biol} and
    \Robject{moltyp} variables Do you notice any problem ?.
  \item Recode the factor levels for the \Robject{mol.biol} and
    \Robject{moltyp} variables using the \Rfunction{factor} function.
  \end{itemize}
\end{Exercise}

\begin{Solution}
<<recodeFactors>>=
levels(psubData$BT)
psubData$BT <-  factor(psubData$BT) 
levels(psubData$BT)
psubData$mol.biol <- factor(psubData$mol.biol)
levels(psubData$mol.biol)
@
\end{Solution}
%$


\section{Compute summary statistics}

\R{} includes several functions that allows you to do a lot while
writing only few lines of code. A good example is the
\Rfunction{aggregate} function that splits the data into subsets,
computes summary statistics for each subset, and returns the result in
a convenient form. For more details regarding this function, please
type in \Rfunction{help("aggregate")} into an \R{} session. We will be
making use of the \Rfunction{formula} and the \Rfunction{data.frame}
methods for the \Rfunction{aggregate} function in this
example. Another useful function for creating contingency tables that
we will be using is the \Rfunction{xtabs} function.

We proceed to create some summary statistics on the variable
\Robject{psubData} using the \Rfunction{aggregate} and
\Rfunction{xtabs} functions.

\begin{Exercise}
  \begin{itemize}
  \item Find the average age of male and females in our subsetted
    metadata variable \Robject{psubData} for the NEG and BCR/ABL
    groups using the \Robject{data.frame} interface for the
    \Rfunction{aggregate} function. The table generated by using the
    \Rfunction{aggregate} should look similar to the one found
    below. Hint: Try passing \Robject{na.rm = TRUE} to the aggregate
    function
    
    \begin{table}[ht]
      \begin{center}
        \begin{tabular}{rllr}
          \hline
          & sex & molBiol & age \\ 
          \hline
          1 & F & BCR/ABL & 39.94 \\ 
          2 & M & BCR/ABL & 40.50 \\ 
          3 & F & NEG & 29.75 \\ 
          4 & M & NEG & 24.86 \\ 
          \hline
        \end{tabular}
      \end{center}
    \end{table}

  \item Recalculate the average age of male and females in our
    subsetted metadata variable \Robject{psubData} for the NEG and
    BCR/ABL groups, this time using the \Robject{formula} interface
    for the \Rfunction{aggregate} function. Make sure that the results
    are identical to the one from the previous step.

  \item The column \Robject{relapse} in \Robject{psubData} is a
    logical vector indicating whether the patient had a relapse or
    not. Create a contingency table of the number of subjects that
    have had a relapse for the samples included in \Robject{psubData}
    using the \Rfunction{xtabs} function and the covariates
    \Robject{relapse}, \Robject{mol.biol} and \Robject{sex}. The table
    generated by using the \Rfunction{xtabs} function should look
    similar to the one found below.

    \begin{table}[ht]
      \begin{center}
        \begin{tabular}{rrr}
          \hline
          & BCR/ABL & NEG \\ 
          \hline
          F & 7.00 & 3.00 \\ 
          M & 9.00 & 18.00 \\ 
          \hline
        \end{tabular}
      \end{center}
    \end{table}

  \end{itemize}
\end{Exercise}

\begin{Solution}
<<contingencytables>>=
aggregate(psubData[, "age", drop = FALSE], 
          by= list(sex= psubData$sex, molBiol= psubData[["mol.biol"]]),  
          FUN = mean, na.rm = TRUE )
aggregate(age ~ sex + mol.biol, data = psubData, FUN = mean)
xtabs(relapse ~ sex + mol.biol, data = psubData)
@
\end{Solution}

\section{Data Visualization}

Base \R{} can produce many different types of statistical
visualizations. Additional packages such as \Rpackage{lattice} or
\Rpackage{ggplot2} extend this functionality. We will explore the
\Rpackage{lattice} package.  A typical call to the \Rpackage{lattice}
function \Rfunction{xyplot} is
<<xyplotFormula, eval=FALSE>>=
xyplot(y ~ x | c, data, groups = g)
@

The arguments to a lattice function can be summarized in terms of
\begin{enumerate}
\item lattice function: A lattice plotting function such as
  \Rfunction{xyplot}, \Rfunction{dotplot} etc.
\item formula: The first argument to a lattice method is a
  formula. The formula for our example is \code{y} \verb+~+ \code{x |
    c}. If the lattice method takes only a single vector as input, the
  formula can be expressed as \verb+~+ \code{ x | c}.
  \begin{itemize}
  \item primary variables: Variables y (Y axis of the plot) and x (X
    axis of the plot) that defines the lattice display separated by
    the \verb+~+ character.
  \item conditioning variable: Variable c in the example separated
    from the primary variables by the character \code{|}. The
    conditioning variable divides the plot into separate panels.
  \end{itemize}
\item grouping variable: The variable g in the example. The grouping
  variable segregates data into subgroups within each panel.
\item data: A \Rclass{data.frame} with column names corresponding to
  the variables y, x, c and g.
\end{enumerate}

\begin{Exercise}

  \begin{itemize}
  \item Load the \Rpackage{lattice} package. Use the
    \Rfunction{bwplot} function to create a box-and-whiskers plot of
    \Robject{age} as a function of \Robject{sex}, conditioning on
    \Robject{mol.biol}.
  \end{itemize}
\end{Exercise}

\begin{Solution}
<<bwplot,fig=TRUE,include=FALSE,prefix=FALSE,width=4,height=4>>=
library(lattice)
plt <- bwplot(age ~ sex | mol.biol, psubData)
print(plt)
@
\end{Solution}

\begin{figure}[p]
  \includegraphics[width=.8\textwidth]{bwplot}
  \caption{Box and whiskers plot summarizing age as a function of sex,
    conditioned on molecular biology.}
\end{figure}


\section{Session information}

<<sessionInfo, results=tex, echo=FALSE>>=
toLatex(sessionInfo())
@

\end{document}

