# BST/STA 226 Discussion 5 February 2007
# Introduction to R

x <- rnorm(30)
y <- rnorm(30)
ls()
x
x[1]
length(x)
mean(x)
help(mean)
sd(x)
class(x)
xy <- cbind(x,y)
xy
class(xy)
dim(xy)
xy[1,1]
xy[1,]
# Get rows with x value greater than 0.5
xy[which(xy[,1]>0.5),]

# Data consists of 26 cases and 500 genes
library(Biobase)
data(sample.exprSet.1)
sample.exprSet.1
class(sample.exprSet.1)
slotNames(sample.exprSet.1)
sample.exprSet.1@exprs
class(expr)
dimnames(expr)
expr <- sample.exprSet.1@exprs
# Expression of one gene for all individuals
expr[1,]
rownames(expr)[1]
hist(expr[1,], main="Expression of First Row", ylab="Frequency", xlab="Expression")

hist(expr[1,1:13], main="Expression of First Row", ylab="Frequency", xlab="Expression", col="blue")
hist(expr[1,14:26], add=T, col="red")

# Expression of all genes for one individuals
expr[,2]
# Get covariate/phenotype data
pData(sample.exprSet.1)
# Get first covariate
sample.exprSet.1@phenoData[[1]]

data(iris)
# Generate means of each attribute
   sapply(1:4, function(x) mean(iris[,x]))

# Generate means of each attribute for each species
   bySpecies <- by(iris[,1:4], iris$Species, mean)
   class(bySpecies)

# Data organized as 50x4x3 array
# Rows are observations, columns are measurements, layers are species
   data(iris3)
   apply(iris3, 2, mean)
   apply(iris3, 3, mean)

   boxplot(iris$Sepal.Length~iris$Species)
   boxplot(iris$Sepal.Length~iris$Species, main="Sepal Length of Irises", xlab="Species", ylab="Sepal Length (mm)")

   plot(Sepal.Length~Petal.Length)
   abline(lm(Sepal.Length~Petal.Length))

   plot(Sepal.Length~Petal.Length, pch="*")
   plot(Sepal.Length~Petal.Length, pch="*", col="blue")
   plot(Sepal.Length~Petal.Length, pch="*", col="blue", cex=2)
   abline(lm(Sepal.Length~Petal.Length), lty=3, col="red", lwd="2")

   split.screen(c(1,2))
   screen(1)
   boxplot(iris$Sepal.Length~iris$Species, main="Sepal Length of Irises", xlab="Species", ylab="Sepal Length (mm)")
   screen(2)
   plot(Sepal.Length~Petal.Length, pch="*")