WIXLIB

xLIST OF .97.107.117.127.137.148.18.28.38.48.58.68.7Mat plot of intensity values for a probe of MLL.B. . . . . . .Density of MLL.B data. . . . . . . . . . . . . . . . . . . . . .Boxplot of the ALL1/AF4 patients. . . . . . . . . . . . . . .Boxplot of the ALL1/AF4 patients after median subtractionand MAD division. . . . . . . . . . . . . . . . . . . . . . . .Venn diagram of seleced ALL genes. . . . . . . . . . . . . . .Boxplot of the ALL1/AF4 patients after median subtractionand MAD division. . . . . . . . . . . . . . . . . . . . . . . . 93. 93. 97. 97. 100. 100Plot of five points to be clustered. . . . . . . . . . . . . . . . . 122Tree of single linkage cluster analysis. . . . . . . . . . . . . . . 122Example of three without clusters. . . . . . . . . . . . . . . . 123Three clusters with different standard deviations. . . . . . . . 123Plot of gene ”CCND3 Cyclin D3” and ”Zyxin” expressions forALL and AML patients. . . . . . . . . . . . . . . . . . . . . . 124Single linkage cluster diagram from gene ”CCND3 Cyclin D3”and ”Zyxin” expressions values. . . . . . . . . . . . . . . . . 124K-means cluster analysis. . . . . . . . . . . . . . . . . . . . . . 126Tree of single linkage cluster analysis. . . . . . . . . . . . . . . 126Plot of kmeans (stars) cluster analysis on CCND3 Cyclin D3and Zyxin discriminating between ALL (red) and AML (black)patients. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 130Vectors of linear combinations. . . . . . . . . . . . . . . . . . . 135First principal component with projections of data. . . . . . . 135Scatter plot of selected genes with row labels on the first twoprincipal components. . . . . . . . . . . . . . . . . . . . . . . 138Single linkage cluster diagram of selected gene expression values.138Biplot of selected genes from the golub data. . . . . . . . . . . 144ROC plot for expression values of CCND3 Cyclin D3. . . . . .ROC plot for expression values of gene Gdf5. . . . . . . . . .Boxplot of expression values of gene a for each leukemia class.Classification tree for gene for three classes of leukemia. . . . .Boxplot of expression values of gene a for each leukemia class.Classification tree of expression values from gene A, B, and Cfor the classification of ALL1, ALL2, and AML patients. . . .Boxplot of expression values from gene CCND3 Cyclin D3 forALL and AML patients . . . . . . . . . . . . . . . . . . . . .149149151151154154156

LIST OF FIGURESClassification tree of expression values from gene CCND3 Cyclin D3 for classification of ALL and AML patients. . . . .8.9 rpart on ALL B-cel 123 data. . . . . . . . . . . . . . . . . .8.10 Variable importance plot on ALL B-cell 123 data. . . . . .8.11 Logit fit to the CCND3 Cyclin D3 expression values. . . . .xi8.89.19.29.3.156159159171G C fraction of sequence ”AF517525.CCND3” along a window of length 50 nt. . . . . . . . . . . . . . . . . . . . . . . . 178Frequency plot of amino acids from accession number AF517525.CCND3.179Frequency plot of amino acids from accession number AL160163.CCND3.17910.1 Graph of probability transition matrix . . . . . . . . . . . . . 19610.2 Evaluation of models by AIC . . . . . . . . . . . . . . . . . . . 21610.3 Tree according to GTR model. . . . . . . . . . . . . . . . . . . 217

xiiLIST OF FIGURES

List of Tables2.1A frequency table and its pie of Zyxin gene. . . . . . . . . . . 183.13.23.3Discrete density and distribution function values of S3 , withp 0.6. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33Built-in-functions for random variables used in this chapter. . 42Density, mean, and variance of distributions used in this chapter. 437.1Data set for principal components analysis.8.18.2Frequencies empirical p-values lower than or equal to 0.01. . . 146Ordered expression values of gene CCND3 Cyclin D3, index2 indicates ALL, 1 indicates AML, cutoff points, number offalse positives, false positive rate, number of true positives,true positive rate. . . . . . . . . . . . . . . . . . . . . . . . . 1709.1BLOSUM50 matrix. . . . . . . . . . . . . . . . . . . . . . . . 186xiii. . . . . . . . . . 134

xivLIST OF TABLES

Chapter 1Brief Introduction into Using RTo get started a gentle introduction to the statistical programming languageR will be given (R Development Core Team, 2009), specific for our purposes.This will solve the practical issues to follow the stream of reasoning. Inparticular, it is briefly explained how to install R and Bioconductor, how toobtain help, and how to perform simple calculations.Since many computations are essentially performed on data vectors, several basic illustrations of this are given. With respect to gene expressions thedata vectors are placed one beneath the other to form a data matrix withthe genes as rows and the patients as columns. The idea of a data matrix isextensively explained and illustrated by several examples. A larger exampleconsists of the classical Golub et al. (1999) data, which will be analyzedfrequently to illustrate statistical procedures.1.1Getting R Started on your PCYou can downloaded R freely from http://cran.r-project.org. Click onyour favorite operating system (Windows, Linux or MacOS) and simply followthe instructions. After a little patience you should be able to start R (Ihaka& Gentleman, 1996) after which a screen is opened with the prompt . Theinput and output of R will be displayed in verbatim typewriting style.All useful functions of R are contained in libraries which are called ”packages”. The standard installation of R makes basic packages available suchas base and stats. From the button Packages at cran.r-project.org itcan be seen that R has a huge number of packages available for a wide scale1

2CHAPTER 1. BRIEF INTRODUCTION INTO USING Rof statistical procedures. To download a specific package you can use thefollowing. install.packages(c("TeachingDemos"),repo "http://cran.r-project.org", dep TRUE)This installs the package TeachingDemos developed by Greg Snow from therepository http://cran.r-project.org. By setting the option dep to TRUEthe packages on which the TeachingDemos depend are also installed. This isstrongly recommended! Alternatively, in the Windows application of R youcan simply click on the Packages button at the top of your screen and followthe instructions. After installing you have to load the package in order to useits functions. For instance, to produce a nice plot of the outcome of throwingtwelve times with a die, you can use the following. library(TeachingDemos) plot(dice(12,1))In the sequel we shall often use packages from Bioconductor, a very usefulopen source software project for the analysis and comprehension of genomicdata. To follow the book it is essential to install Bioconductor on your PCor network. Bioconductor is primarily based on R and can be installed, asfollows. source("http://www.bioconductor.org/biocLite.R") biocLite()Then to download the ALL package from a repository to your system, to loadit, and to make the ALL data (Chiaretti, et. al, 2004) available for usage, youcan use the following. biocLite("ALL") library(ALL) data(ALL)These data will be analyzed extensively later-on in Chapter 5 and 6. Generalhelp on loaded Bioconductor packages becomes available by openVignette().For further information the reader is referred to www.bioconductor.org orto several other URL’s1 .1http://mccammon.ucsd.edu/ bgrant/bio3d/user guide/user cs.conductor

1.2. GETTING HELP3In this and the following chapters we will illustrate many statistical ideasby the Golub et al. (1999) data, see also Section 1.8. The golub data becomeavailable by the following.2 library(multtest) data(golub)R is object-oriented in the sense that everything consists of objects belongingto certain classes. Type class(golub) to obtain the class of the object goluband str(golub) to obtain its structure or content. Type objects() or ls()to view the currently loaded objects, a list probably growing soon to be large.To prevent conflicting definitions, it is wise to remove them all at the end ofa session by rm(list ls()). To quit a session, type q(), or simply click onthe cross in the upper right corner of your screen.1.2Getting helpAll functionalities of R are well-organized in so-called packages. Use the function library() to see which packages are currently installed on your operating system. The packages stats and base are automatically installed, because these contain many basic functionalities. To obtain an overview of thecontent of a package use ls(package:stats) or library(help "stats").Help on the purpose of specific functions can be obtained from the (package)manual by typing a question mark in front of a function. For instance, ?sumgives details on summation. In case you are seeking help on a function whichuses if, simply type apropos("if"). When you are starting with a new concept such as ”boxplot”, it is convenient to have an example showing output(a plot) and programming code. Such is given by example(boxplot). Thefunction history can be useful for collecting previously given commands.Type help.start() to launch an HTML page linking to several wellwritten R manuals such as: ”An Introduction to R”, ”The R Language Definition”, ”R Installation and Administration”, and ”R Data Import/Export”.Further help can be obtained from http://cran.r-project.org. Its ”contributed” page contains well-written freely available on-line books3 and useful reference charts4 . At http://www.r-project.org you can use R site2Functions to read data into R are read.table or read.csv, see also the ”The R DataImport/Export manual”.3”R for Beginners” by Emmanuel Paradis or the ”The R Guide” by Jason Owen4”R reference card” by Tom Short or by Jonathan Baron

4CHAPTER 1. BRIEF INTRODUCTION INTO USING Rsearch, Rseek, or other useful search engines. There are a number of usefulURL’s with information on R.51.3Calculating with RR can be used as a simple calculator. For instance, to add 2 and 3 we simplyinsert the following. 2 3[1] 5In many calculations the natural base e 2.718282 of exponential functionsis used. Such type of functions can be called as follows. exp(1)[1] 2.718282To compute e2 e · e we use exp(2).6 So, indeed, we have ex exp(x), forany value of x.The sum 1 2 3 4 5 can be computed by sum(1:5)[1] 15and the product 5! 5 · 4 · 3 · 2 · 1 by prod(1:5)[1] 1201.4Generating a sequence and a factorIn order to compute so-called quantiles of distributions (see e.g. Section2.1.4) or plots of functions, we need to generate sequences of numbers. Theeasiest way to construct a sequence of numbers is by 1:5[1] 1 2 3 4 55We mention in particular:http://faculty.ucr.edu/ tgirke/Documents/R BioCond/R BioCondManual.html6The argument of functions is always placed between parenthesis ().

1.5. COMPUTING ON A DATA VECTOR5This sequence can also be produced by the function seq, which allows forvarious sizes of steps to be chosen. For instance, in order to compute percentiles of a distribution we may want to generate numbers between zero andone with step size equal to 0.1. seq(0,1,0.1)[1] 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0For plotting and testing of hypotheses we need to generate yet anothertype of sequence, called a “factor”. It is designed to indicate an experimental condition of a measurement or the group to which a patient belongs.7When, for instance, for each of three experimental conditions there are measurements from five patients, the corresponding factor can be generated asfollows. factor - gl(3,5) factor[1] 1 1 1 1 1 2 2 2 2 2 3 3 3 3 3Levels: 1 2 3The three conditions are often called “levels” of a factor. Each of theselevels has five repeats corresponding to the number of observations (patients)within each level (type of disease). We shall further illustrate the idea of afactor soon because it is very useful for purposes of visualization.1.5Computing on a data vectorA data vector is simply a collection of numbers obtained as outcomes frommeasurements. This can be illustrated by a simple example on expressionvalues of a gene. Suppose that gene expression values 1, 1.5, and 1.25 fromthe persons ”Eric”, ”Peter”, and ”Anna” are available. To store these in avector we use the concatenate command c(), as follows. gene1 - c(1.00,1.50,1.25) gene1[1] 1.00 1.50 1.257See e.g. Samuals & Witmer (2003, Chap. 8) for a full explanation of experimentsand statistical principles of design.

6CHAPTER 1. BRIEF INTRODUCTION INTO USING RNow we have created the object gene1 containing three gene expression values. To compute the sum, mean, and standard deviation of the gene expression values we use the corresponding built-in-functions. sum(gene1)[1] 3.75 mean(gene1)[1] 1.25 sum(gene1)/3[1] 1.25 sd(gene1)[1] 0.25 sqrt(sum((gene1-mean(gene1)) 2)/2)[1] 0.25By defining x1 P1.00, x2 1.50, and x3 1.25, the sum of the weightsPcanbe expressed as ni 1 xi 3.75. The mathematical summationisP symbolin R language simply sum. The mean is denoted by x 3i 1 xi /3 1.25and the sample standard deviation asvu 3uXs t (x x)2 /(3 1) 0.25.ii 11.6Constructing a data matrixIn various types of spreadsheets it is custom to store data values in theform of a matrix consisting of rows and columns. In bioinformatics geneexpression values (from several groups of patients) are stored as rows suchthat each row contains the expressions values of the patients correspondingto a particular gene and each column contains all gene expression values fora particular person. To illustrate this by a small example suppose that wehave the following expression values on three genes from Eric, Peter, andAnna.8 gene2 - c(1.35,1.55,1.00) gene3 - c(-1.10,-1.50,-1.25) gene4 - c(-1.20,-1.30,-1.00)8By the function data.entry you can open and edit a screen with the values of amatrix.

1.6. CONSTRUCTING A DATA MATRIX7Before constructing the matrix it is convenient to add the names of the rowsand the columns. To do so we construct the following list. rowcolnames - list(c("gene1","gene2","gene3","gene4"), c("Eric","Peter","Anna"))After the last comma in the first line we give a carriage return for R to comeup with a new line starting with in order to complete a command. Now wecan construct a matrix containing the expression values from our four genes,as follows. gendat - matrix(c(gene1,gene2,gene3,gene4), nrow 4, ncol 3, byrow TRUE, dimnames rowcolnames)Here, nrow indicates the number of rows and ncol the number of columns.The gene vectors are placed in the matrix as rows. The names of the rowsand columns are attached by the dimnames parameter. To see the content ofthe just created object gendat, we print it to the screen. gendatgene1gene2gene3gene4Eric Peter Anna1.00 1.50 1.251.35 1.55 1.30-1.10 -1.50 -1.25-1.20 -1.30 -1.00A matrix such as gendat has two indices [i,j], the first of which refers torows and the second to columns9 . Thus, if you want to print the secondelement of the first row to the screen, then type gendat[1,2]. If you wantto print the first row, then use gendat[1,]. For the second column, usegendat[,2].It may be desirable to write the data to a file for using these in a laterstage or to send these to a college of yours. Consider the following script. write.table(gendat,file "D:/data/gendat.Rdata") gendatread - read.table("D:/data/gendat.Rdata") gendatreadEric Peter Annagene1 1.00 1.50 1.259Indices referring to rows, columns, or elements are always between square brackets [].

8CHAPTER 1. BRIEF INTRODUCTION INTO USING Rgene2 1.35 1.55 1.30gene3 -1.10 -1.50 -1.25gene4 -1.20 -1.30 -1.00An alternative is to use write.csv.101.7Computing on a data matrixMeans or standard deviations of rows or columns are often important fordrawing biologically relevant conclusions. Such type of computations on adata matrix can be accomplished by “for loops”. However, it is much moreconvenient to use the apply functionality on a matrix. To do so we specifythe name of the matrix, indicate rows or columns (1 for rows and 2 forcolumns), and the name of the function. To illustrate this we compute themean of each person (column). apply(gendat,2,mean)Eric PeterAnna0.0125 0.0625 0.0750Similarly, the mean of each gene (row) can be computed. apply(gendat,1,mean)gene1gene2gene3gene41.250000 1.400000 -1.283333 -1.166667It frequently happens that we want to re-order the rows of a matrix accordingto a certain criterion, or, more specifically, the values in a certain columnvector. For instance, to re-order the matrix gendat according to the rowmeans, it is convenient to store these in a vector and to use the functionorder. meanexprsval - apply(gendat,1,mean) o - order(meanexprsval,decreasing TRUE) o[1] 2 1 4 310For more see the ”R Data import/Export” manual, Chapter 3 of the book ”R forBeginners”, or search the internet by the key ”r wiki matrix”.

1.7. COMPUTING ON A DATA MATRIX9Thus gene2 appears first because it has the largest mean 1.4, then gene1with 1.25, followed by gene4 with -1.16 and, finally, gene3 with -1.28. Nowthat we have collected the order numbers in the vector o, we can re-orderthe whole matrix by specifying o as the row index.11 gendat[o,]Eric Peter Annagene2 1.35 1.55 1.30gene1 1.00 1.50 1.25gene4 -1.20 -1.30 -1.00gene3 -1.10 -1.50 -1.25Another frequently occurring problem is that of selecting genes with a certainproperty. We illustrate th

Applied Statistics for Bioinformatics using R Wim P. Krijnen November 10, 2009. ii . all its outcomes are directly available as objects for further programming. R is a rapidly growing language making basic as well as advanced statisti- . to copy and paste scripts from the