# graph a skewed sampling distribution gskew <- function(n = 1, mu = 0, sd = 1, skew = 0) { se <- sd/sqrt(n) x <- seq(mu - 4 * se, mu + 4 * se, length = 501) y <- dskew(x, n, mu, sd, skew) plot(x, y, xlab = "", ylab = "density", main = paste("Distribution of x-bar, sample size =",n), type = "n") lines(x, y, type = "l", col = 4) abline(h = 0) lines(x, dnorm(x, mu, se), col = 8) return(invisible(y)) } # density of skewed distribution dskew <- function(x, n = 1, mu = 0, sd = 1, skew = 0) { se <- sd/sqrt(n) if(skew == 0) y <- dnorm(x, mu, se) else { a <- (skew * sd)/2 b <- mu - (2 * sd)/skew alpha <- 4/skew/skew if(skew > 0) y <- (n * dgamma((n * (x - b))/a, n * alpha))/a else y <- ( - n * dgamma((n * (x - b))/a, n * alpha))/a } return(y) } # density of symmetric bimodal distribution dbimod <- function(x,n=1,mu=0,sd=1,d=1) { if(d >= sd) stop("sd must be larger than d") temp <- matrix(0,nrow=n+1,ncol=length(x)) for(i in 0:n){temp[i+1,] <- dbinom(i,n,.5)*dnorm(x,mu-d+2*i*d/n,sqrt((sd*sd-d*d)/n))} return(apply(temp,2,sum)) } # graph sampling distribution of symmetric bimodal distribution gbimod <- function(n=1, mu=0, sd=1, d=0.5) { if(d >= sd) stop("sd must be larger than d") se <- sd/sqrt(n) x <- seq(mu - 4 * se, mu + 4 * se, length = 201) y <- dbimod(x, n, mu, sd, d) xx <- c(mu - 4 * se, mu + 4 * se) yy <- c(0, max(y, dnorm(mu, mu, se))) plot(xx, yy, xlab = "", ylab = "density", main = paste("Distribution of x-bar, sample size =",n), type = "n") lines(x, y, type = "l", col = 4) abline(h = 0) lines(x, dnorm(x, mu, se), col = 8) return(invisible(y)) } # graph binomial distribution gbinom <- function(n, p, low=0, high=n, scale = F) { sd <- sqrt(n * p * (1 - p)) if(scale && (n > 10)) { low <- max(0, round(n * p - 4 * sd)) high <- min(n, round(n * p + 4 * sd)) } values <- low:high probs <- dbinom(values, n, p) plot(c(low,high), c(0,max(probs)), type = "n", xlab = "Possible Values", ylab = "Probability", main = paste("Binomial Distribution \n", "n =", n, ", p =", p)) lines(values, probs, type = "h", col = 2) abline(h=0,col=3) return(invisible()) }