---
title: "**Bioinformatic Analysis of RNA-sequencing** <br> <br> R Introduction: Files & usefull functions"
subtitle: "Allison Burns, Maxime Jan, Linda Mhalla, Christian Iseli, Nicolas Guex"
output:
  ioslides_presentation:
    css: ../css/styles.css
    incremental: no
    widescreen: true
    toc: yes
    transition: faster
---

# Usefull functions

## Operation on vector {.pre}

:::: {style="display: grid; grid-template-columns: 1fr 1fr; grid-column-gap: 50px; "}
::: {}

Repeat a variable or vector X time

```{r}
MyVector<-c(1,2,3,4,5)
rep(MyVector,times=2)
```


Or you can repeat each element of a vector X time

```{r}
rep(MyVector,each=2)
```


:::
::: {}

You can use paste to add some character (e.g. days)

```{r}
paste("Day",MyVector,sep="_")
```

Or collapse all character into a single character

```{r}
paste("Day",MyVector,sep="_",collapse = ".")
```

:::
::::

## The *apply* functions {.pre}

the apply functions can be usefull to use a function on every item of a vector, list, matrix

:::: {style="display: grid; grid-template-columns: 1fr 1fr; grid-column-gap: 50px; "}
::: {}

- sapply to use a function on item in a vector

```{r}
MyGenes<-c("Arntl","Per2","Cry1","per1")
sapply(MyGenes,toupper)
```

- lapply can be used on a list and will return list

```{r}
lapply(MyGenes,toupper)
```

:::
::: {}

- apply on a matrix

```{r}
MyMatrix<-matrix(c(1,2,3,4,5,6),ncol=2)
MyMatrix

# Compute mean value for each column
apply(MyMatrix,2,mean)
```

To compute for each row, use: apply(MyMatrix,1,mean)

:::
::::

## {.pre}

For sum of matrix, you can simply use RowSum or ColSum as well as rowMeans and colMeans

```{r}
rowSums(MyMatrix)
rowMeans(MyMatrix)
```

```{r}
colSums(MyMatrix)
colMeans(MyMatrix)
```

## *aggregate* function {.pre}

You may want to compute some mean or sd over same samples or genes. You can use the aggregate function

```{r}
MyData<-data.frame(X=1:4,
                   groups=rep(c("A","B"),each=2))
MyData
```


```{r}
aggregate(x = MyData$X, by = list(MyData$groups), FUN = mean) 
```


## Create Functions {.pre}

You can create your own function with some arguments

```{r}
MyFunction<-function(arg1,arg2){
  SumOfArgs<-arg1 + arg2
  return(SumOfArgs)
}
MyFunction(arg1 = 1,arg2 = 10)
```

You can have default arguments

```{r}
MyOtherFunction<-function(arg1=3,arg2="Cortex"){
  RepArg<-rep(arg2,arg1)
  return(RepArg)
}
MyOtherFunction()
```

## {.pre}

You can pass any named or unamed argument using the  '...'

```{r}
MyPlotFunction<-function(x,...){
  
  MyDatasquared<-x^2
  
  plot(MyDatasquared,...)
}
```

The additional arguments are passed to the *plot()* function

:::: {style="display: grid; grid-template-columns: 1fr 1fr; grid-column-gap: 50px; "}
::: {}

```{r,fig.width=4,fig.height=3}
MyPlotFunction(x = rnorm(10) , pch=19, 
               col = "red", main="Dots")
```

:::
::: {}

```{r,fig.width=4,fig.height=3}
MyPlotFunction(x = rnorm(10) , type="l", 
               col="green", main="line")
```


:::
::::


## Local and global variables {.pre}


:::: {style="display: grid; grid-template-columns: 1fr 1fr; grid-column-gap: 50px; "}
::: {}

Within a function, your variable are local. 

I can access within a function a variable defined outside a function (e.g. *Variable_A*)

```{r}
Variable_A<-1


MyFunction<-function(){
  Variable_B<-2
  print(Variable_A)
}

MyFunction()
```

however, Variable_B do not "exist" outside the function

```{r}
exists("Variable_B")
```


:::
::: {}

Also you cannot change Variable_A inside your function with a simple assignment operator '<-'

You need a superassignment operator '<<-'

```{r}
Variable_A<-1
Variable_B<-1
MyFunction<-function(){
  Variable_A<-2
  Variable_B<<-2
  return(invisible(NULL))
}
MyFunction()

Variable_A
Variable_B
```


:::
::::

## Loops {.pre}

:::: {style="display: grid; grid-template-columns: 1fr 1fr; grid-column-gap: 50px; "}
::: {}

You can use for loop

```{r}
for (i in 1:10){
  print(i)
}
```

:::
::: {}

While loop

```{r}
Num<-4
while(Num > 0){
  Num<-Num - 1
  print(Num)
}
```


:::
::::

## if else {.pre}

```{r}
Variable_A<- TRUE
x<-1

if (Variable_A==TRUE){

     x<-x+1

     y<-10

}else{
  x<-x-1
}

x
y
```


## gsub {.pre}

If you want to replace some character into a string, use gsub (see regular expression)

*gsub(pattern, replacement, x)*

```{r}
MyData<-c("Mus Musculus", "Homo Sapiens")
gsub(" ","_",MyData)
```

e.g. a more complex example

```{r}
Pathways<-c("Glycolysis - Rattus norvegicus (rat)")
```


```{r}
gsub(" - Rattus norvegicus \\(rat\\)","",Pathways)
```


## reshape2 package {.pre}

Often, you will recieve data in the form of a data.frame or matrix, and you want to transform it into a data.frame with a single value per row. 

<br> <br>

:::: {style="display: grid; grid-template-columns: 1fr 1fr; grid-column-gap: 50px; "}
::: {}

A usefull package is *reshape2* that will help you to do this:

```{r, echo=F}
MyDataFrame<-data.frame("Genes"=c("Cry1","Per2","Arntl","Cdk4"),
                        "Sample1"=rpois(4,lambda = 4),
                        "Sample2"=rpois(4,lambda = 4),
                        "Sample3"=rpois(4,lambda = 4)
                        )
```


```{r}
MyDataFrame
```

:::
::: {}


We can *melt* this 

```{r}
library(reshape2)
MyNewDF<-melt(MyDataFrame,id.vars = "Genes",
     measure.vars = c("Sample1","Sample2","Sample3"),
     variable.name = "Sample",value.name = "ReadCounts"
     )
MyNewDF
```

:::
::::


# Statistical test

## {.pre}

Many statistical test are available in R: 

- student test: *t.test()*

- non-parametric t-test: *wilcox.test()*

- correlation: *cor.test()* with pearson,spearman or kendall methods

- ANOVA: *aov()*, for type II  or type III ANOVA you can use the the car package (*car::Anova()*)

- ANOVA post-hoc tukey using *TukeyHSD()*

- Linear model using *lm()* (we will see how to write your model !)

- Mixed effect model using lme4 package


## Linear Regression {.pre}

:::: {style="display: grid; grid-template-columns: 1fr 1fr; grid-column-gap: 50px; "}
::: {}


Linear model is written in this form: 

$y = \beta_1X_1 + \beta_2X_2 + \beta_3X_3 ...$

- where y is your variable of interest (dependent variable)

- X your predictors (independent variables)

- $\beta$ values are the weights of your predictors (what you generally want to test)

```{r,echo=F}
set.seed(42)
Age<-rpois(20,20)
Gene_Expression<-20+2*Age+rnorm(20,sd=4)
```



```{r,fig.width=3,fig.height=3}
plot(Age,Gene_Expression,pch=19)
```


:::
::: {}

We test for an effect of age on gene expression

```{r}
 summary(lm(Gene_Expression~Age))[[4]]
```


And found a significant effect of age: 

$GeneExpression = \beta_0 + \beta_1 * Age + \epsilon$; $\beta_1=$ `r coef(lm(Gene_Expression~Age))[2]``
```{r,fig.width=3,fig.height=3,echo=F}
plot(Age,Gene_Expression,pch=19)
abline(lm(Gene_Expression~Age),
       col="darkblue",
       lwd=2)
```


:::
::::

## Design Matrix {.pre}

Depending on your design, DEseq2 may require a design matrix and which coefficient you want to test.


:::: {style="display: grid; grid-template-columns: 1fr 1fr; grid-column-gap: 50px; "}
::: {}

Your design may look like this:

```{r,echo=F}
MyMetadata<-data.frame("Gene_Expression"=rnorm(9),
                       "Treatment"=rep(c("Control","A_vehic","Treatment"),each=3),
                       "Batch" = rep(c("B1","B2","B3"),3))
MyMetadata
```

:::
::: {}

<br>

You can write a model matrix that will look like this, but your Treatment *A_vehic* becomes the reference (your coefficient 4 represent the difference of $Control - A\_vehic$, and coefficient 5: $Treatment - A\_vehic$)

```{r}
head(model.matrix(~Batch + Treatment,data=MyMetadata),9)
```

:::
::::


## {.pre}

To solve this problem, write Treatment as factor and relevel it:

```{r}
MyMetadata$Treatment<-factor(MyMetadata$Treatment,
                             levels = c("Control","A_vehic","Treatment"))
```


Now your coefficient 4 represent the difference of $A\_vehic - Contro l$, and coefficient 5: $Treatment - Control$)

```{r}
head(model.matrix(~Batch + Treatment,data=MyMetadata),9)
```

# Files

## Read files

- You can read many different files in R (csv, txt, xlsx etc...)

- If you file has a header, add *header=T*

- You can use the built-in function like *read.table()* or *read.csv()* for text files

- The package openxlsx will read excel format files

- You can also use the *Import Dataset* in Rstudio (upper-right corner)

- You can press *tab* while writing the path to a file to help you

- You can also skip some line in a file using the *skip* option in *read.table()*

## Write files {.pre}

- To write files use *write.table()* function.

- Some option are very usefull, like:

```{r eval=F}
write.table(MyDataFrame,
            quote=F, # Avoid adding "" to character
            sep="\t", # Use tab to separate data
            col.names = T,# Add colnames to the file
            row.names = T # Add rownames to the file
            )
```

or save directly your variable in a Rdata object

```{r,eval=F}
save(MyDataFrame,MyDataFrame2,file="MyDataframes.Rdata")
```

and load the file using *load()*

```{r, eval=F}
load("MyDataframes.Rdata")
```





