---
title: "Contrôle négatif pour test de Welch"
author: "Jacques van Helden"
date: "11/6/2017"
output:
  html_document: 
    fig_caption: yes
    highlight: zenburn
    theme: cerulean
    toc: yes
    toc_depth: 3
    toc_float: yes
    code_folding: hide
  pdf_document: default
  word_document: default
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```

## The Normal distribution


### Normal density

```{r}
# ?Normal
x <- seq(from=-5, to=5, by=0.1)
plot(x, y = dnorm(x = x, mean = 0, sd = 1), type="l",
     main="Normal density", xlab="X", ylab="dnorm(0, 1, x)", col="blue")
grid()
abline(v=0)
```

### Normal probability

The function `pnorm(q=x, m, s)` gives the probability to observe a value smaller than or equal to $x$: $P(X \le x)$. 

```{r}
plot(x, y = pnorm(q = x, mean = 0, sd = 1), type="l",
     main="Normal probability (CDF)", xlab="X", ylab="dnorm(0, 1, x)", col="blue")
grid()
abline(v=0)
abline(h=0.5, col="darkgreen")

```

The function `pnorm(q=x, m, s, lower.tail=FALSE)` gives the probability to observe a value greater than to $x$: $P(X \gt x)$. 
```{r}
plot(x, y = pnorm(q = x, mean = 0, sd = 1, lower.tail = FALSE), type="l",
     main="Normal probability (upper tail)", xlab="X", ylab="dnorm(0, 1, x)", col="blue")
grid()
abline(v=0)
abline(h=0.5, col="darkgreen")

```

### qnorm

```{r}
qnorm(p = 0.05, mean = 0, sd = 1)
qnorm(p = 0.95, mean = 0, sd = 1)
```


### rnorm

Generates random numbers following a normal distribution.

```{r}
rnorm(n = 10, mean = 0, sd = 1)
```

We can draw an hitogram of $n=10,000$ values in order to check if it has a normal shape. 


```{r}
hist(rnorm(n=10000, mean = 0, sd = 1), breaks=50)
```



## Generating a negative congtrol set

We will generate two data frames to simulate a microarray data set with $g=10000$ genes with respectively $n_1 = n_2 = 40$ samples for groups 1 and 2,  where each gene (row) is under null hypothesis $H_0: \mu_1=\mu_2$. 

```{r parameters}
g <- 10000

## First group
n1 <- 40
mu1 <- 0
sigma1 <- 1

## Second group
n2 <- 40
mu2 <- 0
sigma2 <- 2
  
```

```{r generate_random_data}
group1 <- matrix(nrow=g, rnorm(n = g * n1, mean = mu1, sd = sigma1))
# View(group1)

group2 <- matrix(nrow=g, rnorm(n = g * n2, mean = mu2, sd = sigma2))
dim(group2)

## Merge the two datasets in a single matrix
rnorm.data <- cbind(group1, group2)
class(rnorm.data)
dim(rnorm.data)

```


```{r mean_variance_plots}
## Mean-variance plot per group
plot(apply(group2, 1, mean), 
     apply(group2, 1, sd), col="purple", 
     ylim=c(0, 4))
lines(apply(group1, 1, mean), 
     apply(group1, 1, sd), col="darkgreen", type="p")
grid()

abline(v=mu1)
abline(h=sigma1)
abline(h=sigma2)
```

### Mean comparison plot

```{r mean_comparison, fig.cap="Comparison of sample means between the two groups. "}
plot(apply(group1, 1, mean),
     apply(group2, 1, mean),
     col="grey"
     )

```



## Apply Welch test to each row of the simulated dataset

## Interpretation of the results of the negative controls


## Conclusions and perspectives


## References


