---
title: "Lab 4- Transform a discrete variable and generate from cdf"
author: "Amber Hackstadt"
date: "`r Sys.Date()`"
output: html_document
---

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

# Generate a continuous random variable

###  Drawing from exponential distribution with $\lambda = 1$

#### Use inverse of cdf
```{r}
n<-1000
lambda = 1
y = runif(n) ## generate n numbers between 0 and 1
f.inv = function(x){-1/lambda*log(1-x)}
X.exp = f.inv(y)
```

#### Compare empirical density to true density using the dexp function in R
```{r}
 ## empirical density
hist(X.exp,freq=FALSE, main  = "Inverse cdf density vs Truth")    ## empirical density

## true density
t<-c(1:100)/10 
lines(t,dexp(t,1),col = "red") 
```

#### Compare the empirical cdf to the true cdf
```{r}
plot(ecdf(X.exp),main = "rexp, Empirical cdf vs truth")                  ## empirical cdf using ecdf function
lines(t,pexp(t,1),col=2)    ## true cdf using pexp function
```

#### Compare to using the rexp function in R
```{r}
X.exp2 = rexp(n,1)
hist(X.exp2,freq=FALSE, main = "rexp vs truth")
t<-c(1:100)/10 
lines(t,dexp(t,1),col = "red") 

plot(ecdf(X.exp2),main = "rexp, Empirical cdf vs truth")                  ## empirical cdf versus true cdf
lines(t,pexp(t,1),col=2)    ## true cdf using pexp function
```

#### Check probability-integral transformation: Drawing X from Exp(1)), computing Y=F(X), and show that Y~UNIF(0,1)
```{r}
y<-pexp(X.exp,1) ### compute F(x)
hist(y,freq=FALSE)  ## look at pdf of Y
s<-c(0:100)/100
lines(s,dunif(s,0,1),col = "red") ### compare to uniform (0,1) 
```

###  Drawing from a chi-square distribution with 1 degree of freedom

#### Generate 10000 from chi-square(1) distribution. Compare empirical density and cdf with true density and cdf. 
```{r}
n<-10000
X<-rchisq(n,1)
hist(X,freq=FALSE, main = "Chi-squared empirical cdf with true density")
t<-c(0:1200)/100
lines(t,dchisq(t,1),col = 2)  ## empirical cdf with true density

plot(ecdf(X), main = "Chi-squared empirical cdf with true cdf")
lines(t,pchisq(t,1),col=2)  ## empirical cdf with true cdf
```


#### Generate 10000 U from UNIF(0,1) and then generate V from chi-squared using the inverse cdf of a chi-squared distribution.
```{r}
n <- 1000
U <- runif(n)
V = qchisq(U,1)
hist(V,  freq = FALSE,main = "Chi-squared using inverse cdf, empirical cdf with true density")
t<-c(0:2000)/100
lines(t,dchisq(t,1), col = "red")  ## empirical pdf with true density

plot(ecdf(X), main = "Chi-squared using inverse cdf, empirical cdf with true cdf")
t<-c(0:2000)/100
lines(t,pchisq(t,1),col=2)  ## empirical cdf with true cdf
```


```{r}
hist(X,main = "compare both runif and inverse cdf",freq = FALSE, col = "red")
hist(V, col = "grey",add =TRUE, freq = FALSE)
```

# Simulation to estimate the pmf for a transformed discrete random variable
X  is the value shown when you roll a six sided die and $Y = (X-2)^2$. Find f_Y(y).

### Simulation code
```{r }
## ++++++++++++++++
## simulation parameters
## ++++++++++++++++
num.sims<-10000
seed.use <- 124

set.seed(seed.use)

## ++++++++++++++++
## create place holders
## ++++++++++++++++
Y.vals <- NULL

## ++++++++++++++++
## run simulation
## ++++++++++++++++
for(i in 1:num.sims) {
  x<-sample(1:6,1,replace=TRUE)
  Y.vals[i] = (x-2)^2
}
```

### Summarize results
```{r }
## ++++++++++++++++
## calculate probability of all different birthday
## ++++++++++++++++
## simulation result
sum.results = data.frame(round(table(Y.vals)/num.sims,3))
colnames(sum.results) = c("Y","P(Y=y)")
rownames(sum.results) =NULL

pander(sum.results)
```