} else if (type == "prop" || type == "proportion")
{
for (i in 1:nReps){
ci[i, 1] <- mean(x) - qnorm(0.975)*sqrt((mean(x)*(1-mean(x)))/10000)
ci[i, 2] <- mean(x) + qnorm(0.975)*(sd(x)/sqrt(nReps))
cover[i] <- ifelse(ci[i,1] < 0.7 && ci[i,2] > 0.7, 1, 0)
}
}
mean(cov)
}
coverage(nReps = 10, type = "mean")
coverage <- function(nReps = 10000, type) #where x is the vector containing the data
# (nReps in length)
{
ci <- matrix(nrow = nReps, ncol = 2)
cover <- NULL
if (type == "mean")
{
for(i in 1:nReps){
x <- avg.kids()
ci[i, 1] <- mean(x) - qnorm(0.975)*(sd(x)/sqrt(nReps))
ci[i, 2] <- mean(x) + qnorm(0.975)*(sd(x)/sqrt(nReps))
cover[i] <- ifelse((ci[i,1] < 75 && ci[i,2]) > 75, 1, 0)
}
} else if (type == "prop" || type == "proportion")
{
for (i in 1:nReps){
ci[i, 1] <- mean(x) - qnorm(0.975)*sqrt((mean(x)*(1-mean(x)))/10000)
ci[i, 2] <- mean(x) + qnorm(0.975)*(sd(x)/sqrt(nReps))
cover[i] <- ifelse(ci[i,1] < 0.7 && ci[i,2] > 0.7, 1, 0)
}
}
mean(cover)
}
coverage(nReps = 10, type = "mean")
ci <- matrix(nrow = nReps, ncol = 2)
cover <- NULL
for(i in 1:nReps){
x <- avg.kids()
ci[i, 1] <- mean(x) - qnorm(0.975)*(sd(x)/sqrt(nReps))
ci[i, 2] <- mean(x) + qnorm(0.975)*(sd(x)/sqrt(nReps))
cover[i] <- ifelse((ci[i,1] < 75 && ci[i,2]) > 75, 1, 0)
}
nReps <- 10
for(i in 1:nReps){
x <- avg.kids()
ci[i, 1] <- mean(x) - qnorm(0.975)*(sd(x)/sqrt(nReps))
ci[i, 2] <- mean(x) + qnorm(0.975)*(sd(x)/sqrt(nReps))
cover[i] <- ifelse((ci[i,1] < 75 && ci[i,2]) > 75, 1, 0)
}
cover <- NULL
for(i in 1:nReps){
x <- avg.kids()
ci[i, 1] <- mean(x) - qnorm(0.975)*(sd(x)/sqrt(nReps))
ci[i, 2] <- mean(x) + qnorm(0.975)*(sd(x)/sqrt(nReps))
cover[i] <- ifelse((ci[i,1] < 75 && ci[i,2]) > 75, 1, 0)
}
ci <- matrix(nrow = nReps, ncol = 2)
for(i in 1:nReps){
x <- avg.kids()
ci[i, 1] <- mean(x) - qnorm(0.975)*(sd(x)/sqrt(nReps))
ci[i, 2] <- mean(x) + qnorm(0.975)*(sd(x)/sqrt(nReps))
cover[i] <- ifelse((ci[i,1] < 75 && ci[i,2]) > 75, 1, 0)
}
cover
for(i in 1:nReps){
x <- avg.kids()
ci[i, 1] <- mean(x) - qnorm(0.975)*(sd(x)/sqrt(nReps))
ci[i, 2] <- mean(x) + qnorm(0.975)*(sd(x)/sqrt(nReps))
cover[i] <- ifelse((ci[i,1] < 75 && ci[i,2] > 75), 1, 0)
}
cover
rm(ci, cover, i , nReps, x)
coverage(type = "prop")
x <- prop.kids()
ci[i, 1] <- mean(x) - qnorm(0.975)*sqrt((mean(x)*(1-mean(x)))/10000)
x <- prop.kids()
rm9x
rm(x)
coverage(type = "prop")
coverage <- function(nReps = 10000, type) #where x is the vector containing the data
# (nReps in length)
{
ci <- matrix(nrow = nReps, ncol = 2)
cover <- NULL
if (type == "mean")
{
for(i in 1:nReps){
x <- avg.kids()
ci[i, 1] <- mean(x) - qnorm(0.975)*(sd(x)/sqrt(nReps))
ci[i, 2] <- mean(x) + qnorm(0.975)*(sd(x)/sqrt(nReps))
cover[i] <- ifelse((ci[i,1] < 75 && ci[i,2] > 75), 1, 0)
}
} else if (type == "prop" || type == "proportion")
{
for (i in 1:nReps){
x <- prop.kids()
ci[i, 1] <- mean(x) - qnorm(0.975)*sqrt((mean(x)*(1-mean(x)))/10000)
ci[i, 2] <- mean(x) + qnorm(0.975)*(sd(x)/sqrt(nReps))
cover[i] <- ifelse(ci[i,1] < 0.7 && ci[i,2] > 0.7, 1, 0)
}
}
mean(cover)
}
coverage(type = "prop")
avg.kids.ci <- function(nReps = 10000)
{
tot_family <- NULL #holds how many kids each family has, including the boy. e.g. if the family has
# 4 children and the 4th is the boy, they had 4 children to get to the boy; will contain nReps amount
# of values
for (i in 1:nReps)
{
prob.boy <- runif(1,0,1) #pulls the probability for the i-th family to have a boy
prob.girl <- 1- prob.boy
count <- 0
family <- NULL #holds what kids get born in each family
while(sum(family) == 0)
{
count <- count + 1
born <- sample(0:1, 1, prob = c(prob.girl, prob.boy)) # girl = 0, boy = 1
family[count] <- born
}
tot_family[i] <- length(family)
}
c(mean(tot_family) + c(-1,1)*qnorm(0.975)*(sd(tot_family)/sqrt(nReps)))
}
avg.kids.ci()
prop.kids.ci <- function(nReps = 10000)
{
tot_family <- NULL
for (i in 1:nReps)
{
prob.boy <- runif(1,0,1) #pulls the probability for the i-th family to have a boy
prob.girl <- 1- prob.boy
born <- sample(0:1, 2, replace = TRUE, prob = c(prob.girl, prob.boy))
tot_family[i] <- ifelse(sum(born) >= 1, 1, 0)
}
c(mean(tot_family) + c(-1,1)*qnorm(0.975)*sqrt((mean(tot_family*(1-mean(tot_family)))/10000))
}
c(mean(tot_family) + c(-1,1)*qnorm(0.975)*sqrt((mean(tot_family*(1-mean(tot_family)))/10000))
}
prop.kids.ci
# find the coverage of both of these situations if the true mean is 75 and the true proportion is 0.7
# and assess monte carlo error
#***** arbitrarily chose 0.7
coverage <- function(nReps = 10000, type) #where x is the vector containing the data
# (nReps in length)
{
ci <- matrix(nrow = nReps, ncol = 2)
cover <- NULL
if (type == "mean")
{
for(i in 1:nReps){
x <- avg.kids()
ci[i, 1] <- mean(x) - qnorm(0.975)*(sd(x)/sqrt(nReps))
ci[i, 2] <- mean(x) + qnorm(0.975)*(sd(x)/sqrt(nReps))
cover[i] <- ifelse((ci[i,1] < 75 && ci[i,2] > 75), 1, 0)
}
} else if (type == "prop" || type == "proportion")
{
for (i in 1:nReps){
x <- prop.kids()
ci[i, 1] <- mean(x) - qnorm(0.975)*sqrt((mean(x)*(1-mean(x)))/10000)
ci[i, 2] <- mean(x) + qnorm(0.975)*(sd(x)/sqrt(nReps))
cover[i] <- ifelse(ci[i,1] < 0.7 && ci[i,2] > 0.7, 1, 0)
}
}
mean(cover)
}
prop.kids.ci
prop.kids.ci()
prop.kids.ci <- function(nReps = 10000)
{
tot_family <- NULL
for (i in 1:nReps)
{
prob.boy <- runif(1,0,1) #pulls the probability for the i-th family to have a boy
prob.girl <- 1- prob.boy
born <- sample(0:1, 2, replace = TRUE, prob = c(prob.girl, prob.boy))
tot_family[i] <- ifelse(sum(born) >= 1, 1, 0)
}
c(mean(tot_family) + c(-1,1)*qnorm(0.975)*sqrt((mean(tot_family*(1-mean(tot_family)))/10000))
}
prop.kids.ci <- function(nReps = 10000)
{
tot_family <- NULL
for (i in 1:nReps)
{
prob.boy <- runif(1,0,1) #pulls the probability for the i-th family to have a boy
prob.girl <- 1- prob.boy
born <- sample(0:1, 2, replace = TRUE, prob = c(prob.girl, prob.boy))
tot_family[i] <- ifelse(sum(born) >= 1, 1, 0)
}
c(mean(tot_family) + c(-1,1)*qnorm(0.975)*sqrt((mean(tot_family*(1-mean(tot_family)))/10000))
}
c(mean(tot_family) + c(-1,1)*qnorm(0.975)*sqrt((mean(tot_family*(1-mean(tot_family)))/10000)))
prop.kids.ci <- function(nReps = 10000)
{
tot_family <- NULL
for (i in 1:nReps)
{
prob.boy <- runif(1,0,1) #pulls the probability for the i-th family to have a boy
prob.girl <- 1- prob.boy
born <- sample(0:1, 2, replace = TRUE, prob = c(prob.girl, prob.boy))
tot_family[i] <- ifelse(sum(born) >= 1, 1, 0)
}
c(mean(tot_family) + c(-1,1)*qnorm(0.975)*sqrt((mean(tot_family*(1-mean(tot_family)))/10000)))
}
prop.kids.ci
prop.kids.ci()
coverage <- function(nReps = 10000, type) #where x is the vector containing the data
# (nReps in length)
{
ci <- matrix(nrow = nReps, ncol = 2)
cover <- NULL
if (type == "mean")
{
for(i in 1:nReps){
x <- avg.kids()
cover[i] <- ifelse((x[1] < 75 && x[2] > 75), 1, 0)
}
} else if (type == "prop" || type == "proportion")
{
for (i in 1:nReps){
x <- prop.kids()
cover[i] <- ifelse(x[1] < 0.7 && x[2] > 0.7, 1, 0)
}
}
mean(cover)
}
coverage(type = "mean")
coverage <- function(nReps = 10000, type) #where x is the vector containing the data
# (nReps in length)
{
ci <- matrix(nrow = nReps, ncol = 2)
cover <- NULL
if (type == "mean")
{
for(i in 1:nReps){
x <- avg.kids.ci()
cover[i] <- ifelse((x[1] < 75 && x[2] > 75), 1, 0)
}
} else if (type == "prop" || type == "proportion")
{
for (i in 1:nReps){
x <- prop.kids.ci()
cover[i] <- ifelse(x[1] < 0.7 && x[2] > 0.7, 1, 0)
}
}
mean(cover)
}
coverage(type = "mean")
help(sapply)
sapply(1:10000, avg.kids.ci)
#on a certain planet, all the families want a boy. they will keep having children until they get a boy.
# what is the mean number of children it takes to get a boy (including the boy)
# if the chance of having a boy in each family follows a uniform distribution between 0 and 1?
# In essence, the probability of having a boy for each family is different and follows the previously
# stated distribution. don't forget to assess monte carlo error.
#******* give how to generate random numbers from a uniform distribution? we can also use a different
# distribution
nReps <- 1000
avg.kids.ci <- function(nReps = 1000)
{
tot_family <- NULL #holds how many kids each family has, including the boy. e.g. if the family has
# 4 children and the 4th is the boy, they had 4 children to get to the boy; will contain nReps amount
# of values
for (i in 1:nReps)
{
prob.boy <- runif(1,0,1) #pulls the probability for the i-th family to have a boy
prob.girl <- 1- prob.boy
count <- 0
family <- NULL #holds what kids get born in each family
while(sum(family) == 0)
{
count <- count + 1
born <- sample(0:1, 1, prob = c(prob.girl, prob.boy)) # girl = 0, boy = 1
family[count] <- born
}
tot_family[i] <- length(family)
}
c(mean(tot_family) + c(-1,1)*qnorm(0.975)*(sd(tot_family)/sqrt(nReps)))
}
avg.kids.ci()
# on the same planet on a different island, all the families also want a boy child.
# however, due to the size of the island, all families can only have a max of two children. What
# proportion of families actually have a boy child? the same probability situation from the above
# problem still applies. don't forget to assess monte carlo error.
prop.kids.ci <- function(nReps = 1000)
{
tot_family <- NULL
for (i in 1:nReps)
{
prob.boy <- runif(1,0,1) #pulls the probability for the i-th family to have a boy
prob.girl <- 1- prob.boy
born <- sample(0:1, 2, replace = TRUE, prob = c(prob.girl, prob.boy))
tot_family[i] <- ifelse(sum(born) >= 1, 1, 0)
}
c(mean(tot_family) + c(-1,1)*qnorm(0.975)*sqrt((mean(tot_family*(1-mean(tot_family)))/10000)))
}
prop.kids.ci()
# find the coverage of both of these situations if the true mean is 75 and the true proportion is 0.7
# and assess monte carlo error
#***** arbitrarily chose 0.7
coverage <- function(nReps = 1000, type) #where x is the vector containing the data
# (nReps in length)
{
ci <- matrix(nrow = nReps, ncol = 2)
cover <- NULL
if (type == "mean")
{
for(i in 1:nReps){
x <- avg.kids.ci()
cover[i] <- ifelse((x[1] < 75 && x[2] > 75), 1, 0)
}
} else if (type == "prop" || type == "proportion")
{
for (i in 1:nReps){
x <- prop.kids.ci()
cover[i] <- ifelse(x[1] < 0.7 && x[2] > 0.7, 1, 0)
}
}
mean(cover)
}
coverage(type = "mean")
#on a certain planet, all the families want a boy. they will keep having children until they get a boy.
# what is the mean number of children it takes to get a boy (including the boy)
# if the chance of having a boy in each family follows a uniform distribution between 0 and 1?
# In essence, the probability of having a boy for each family is different and follows the previously
# stated distribution. don't forget to assess monte carlo error.
#******* give how to generate random numbers from a uniform distribution? we can also use a different
# distribution
nReps <- 1000
avg.kids.ci <- function(nReps = 1000)
{
tot_family <- NULL #holds how many kids each family has, including the boy. e.g. if the family has
# 4 children and the 4th is the boy, they had 4 children to get to the boy; will contain nReps amount
# of values
for (i in 1:nReps)
{
prob.boy <- runif(1,0,1) #pulls the probability for the i-th family to have a boy
prob.girl <- 1- prob.boy
count <- 0
family <- NULL #holds what kids get born in each family
while(sum(family) == 0)
{
count <- count + 1
born <- sample(0:1, 1, prob = c(prob.girl, prob.boy)) # girl = 0, boy = 1
family[count] <- born
}
tot_family[i] <- length(family)
}
c(mean(tot_family) + c(-1,1)*qnorm(0.975)*(sd(tot_family)/sqrt(nReps)))
}
avg.kids.ci()
# on the same planet on a different island, all the families also want a boy child.
# however, due to the size of the island, all families can only have a max of two children. What
# proportion of families actually have a boy child? the same probability situation from the above
# problem still applies. don't forget to assess monte carlo error.
prop.kids.ci <- function(nReps = 1000)
{
tot_family <- NULL
for (i in 1:nReps)
{
prob.boy <- runif(1,0,1) #pulls the probability for the i-th family to have a boy
prob.girl <- 1- prob.boy
born <- sample(0:1, 2, replace = TRUE, prob = c(prob.girl, prob.boy))
tot_family[i] <- ifelse(sum(born) >= 1, 1, 0)
}
c(mean(tot_family) + c(-1,1)*qnorm(0.975)*sqrt((mean(tot_family*(1-mean(tot_family)))/10000)))
}
prop.kids.ci()
# find the coverage of both of these situations if the true mean is 75 and the true proportion is 0.7
# and assess monte carlo error
#***** arbitrarily chose 0.7
coverage <- function(nReps = 1000, type) #where x is the vector containing the data
# (nReps in length)
{
ci <- matrix(nrow = nReps, ncol = 2)
cover <- NULL
if (type == "mean")
{
for(i in 1:nReps){
x <- avg.kids.ci()
cover[i] <- ifelse((x[1] < 75 && x[2] > 75), 1, 0)
}
} else if (type == "prop" || type == "proportion")
{
for (i in 1:nReps){
x <- prop.kids.ci()
cover[i] <- ifelse(x[1] < 0.7 && x[2] > 0.7, 1, 0)
}
}
mean(cover)
}
coverage(type = "mean")
x <- avg.kids.ci(10000000)
x <- avg.kids.ci(1000000)
x
(x[2]-x[1])2
(x[2]-x[1])/2
(x[2]-x[1])
(x[2]-x[1])/2
((x[2]-x[1])/2) + x[1]
13.81109 - ((x[2]-x[1])/2)
13.81109 + ((x[2]-x[1])/2)
x <- avg.kids.ci(1000000)
x <- avg.kids.ci(1000000)
13.81104
x
(x[2]-x[1])/2
(x[2]-x[1])/2 + x[1]
(x[2]-x[1])/2 - x[2]
data(galaxies)
xbar <- zn <- rep(NA, nsim)
#CLT simulation
nsim <- 100000
#exponential(5)
mu <- 1/5
sig2 <- (1/5)^2
sig <- sqrt(sig2)
xbar <- zn <- rep(NA, nsim)
n <- 5
for(i in 1:nsim)
{
x <- rexp(n, 5)
xbar[i] <- mean(x)
zn[i] <- (mean(x)-mu)/(sqrt(sig2/n))
}
plot(density(zn)) #smooth histogram
#distribution of zn
curve(dnorm(x, 0, 1))
#distribution of zn
curve(dnorm(x, 0, 1), -3, 3, add = TRUE, col = 'red')
plot(density(zn)) #smooth histogram/density
#distribution of zn
curve(dnorm(x, 0, 1), -3, 3, add = TRUE, col = 'red')
#bernoulli
mu <- .1
sig2 <- .1*.9
sig <- sqrt(sig2)
xbar <- zn <- rep(NA, nsim)
n <- 100 #need at least a hundred to make it look good.
for(i in 1:nsim)
{
x <- rbinom(n, 1, .1)
xbar[i] <- mean(x)
zn[i] <- (mean(x)-mu)/(sqrt(sig2/n))
}
plot(density(zn)) #smooth histogram/density
#distribution of zn; sampling distribution
curve(dnorm(x, 0, 1), -3, 3, add = TRUE, col = 'red')
#standard normal distribution
load("C://Users//angel//Documents//PPM")
load("C://Users//angel//Documents//PPM//PPM.Rproj")
load("C:\\Users\\angel\\Documents\\PPM\\PPM.Rproj"")
""
load("C:\\Users\\angel\\Documents\\PPM\\PPM.Rproj")
setwd("C:\\Users\\angel\\Documents\\PPM\\R")
load("PPM.Rproj")
load(PPM.Rproj)
load("PPM.Rproj")
library(PPM, lib.loc = "C:\\Users\\angel\\Documents\\PPM")
library(PPM, lib.loc = "C:\\Users\\angel\\Documents\\PPM\\")
library("PPM", lib.loc = "C:\\Users\\angel\\Documents\\PPM\\")
library("PPM")
help(install.packages)
install.packages("C:\\Users\\angel\\Documents\\PPM\\PPM.Rproj)
""
"
install.packages("C:\\Users\\angel\\Documents\\PPM\\PPM.Rproj")
library(PPM)
library(PPM.Rproj)
library("PPM.Rproj")
install.packages("C:\\Users\\angel\\Documents\\PPM\\PPM)
"
)
install.packages("C:\\Users\\angel\\Documents\\PPM\\PPM")
library(PPM)
installed.packages
installed.packages()
install.packages("C:\\Users\\angel\\Documents\\PPM\\PPM.Rproj)
"""))
install.packages("C:\\Users\\angel\\Documents\\PPM\\PPM.Rproj")
install.packages("C:\\Users\\angel\\Documents\\PPM\\PPM")
install("C:\\Users\\angel\\Documents\\PPM\\PPM")
"install.packages("C://Users//angel//Documents//PPM//PPM.Rproj",repos=NULL, type="source",repos=NULL, type="source")
install.packages("C://Users//angel//Documents//PPM//PPM.Rproj",repos=NULL, type="source",repos=NULL, type="source")
install.packages("C://Users//angel//Documents//PPM//PPM.Rproj", type="source",repos=NULL, type="source")
install.packages("C://Users//angel//Documents//PPM//PPM.Rproj", repos=NULL, type="source")
installed.packages()
load("C://Users//angel//Documents//PPM//PPM.Rproj")
install.packages("C://Users//angel//Documents//PPM//PPM.Rproj",repos=NULL, type="source")
install.packages("C://Users//angel//Documents//PPM.zip",repos=NULL, type="source")
installed.packages()
rppm(0,1)
library("PPM")
install.packages("C://Users//angel//Documents//PPM.zip",repos=NULL, type="source")
load("C://Users//angel//Documents//PPM.zip")
install.packages("C:\\Users\\angel\\Documents\\PPM.zip")
install.packages("~/PPM.zip", repos = NULL, type = "win.binary")
library(PPM)
library("PPM")
install.packages("~/PPM.zip", repos = NULL, type = "win.binary")
library("PPM")
installed.packages()
install.packages.zip()
help(install.packages.zip)
??install.packages.zip
install.packages("MCMCpack"
)
library("MCMCpack")