In ddarmon/confdist: confdist, a package for creating confidence functions

weightchange <- read.csv("weight_change_data.csv")

Low Carb

Conf Dist

xbar.lowcarb <- mean(weightchange$weightchange[weightchange$diet=="Low Carb"])
s.lowcarb <- sd(weightchange$weightchange[weightchange$diet=="Low Carb"])
n.lowcarb <- length(weightchange$weightchange[weightchange$diet=="Low Carb"])
xbar.lowcarb
s.lowcarb
n.lowcarb

lowcarb.dist <- curve(1-pt((xbar.lowcarb - mu)/(s.lowcarb/sqrt(n.lowcarb)), df = (n.lowcarb-1)), from = -20, to = 0, xname = 'mu')

Checking Conf Dist

?t.test
t.test(weightchange$weightchange[weightchange$diet=="Low Carb"], alternative = "greater", mu=xbar.lowcarb)

curve(1-pt((xbar.lowcarb - mu)/(s.lowcarb/sqrt(n.lowcarb)), df = (n.lowcarb-1)), from = -20, to = 0, xname = 'mu')
  low.carb.t.test <- t.test(weightchange$weightchange[weightchange$diet=="Low Carb"], alternative = "greater", mu=-12.5)


abline(h=low.carb.t.test$p.value, v=-12.5)

curve(1-pt((xbar.lowcarb - mu)/(s.lowcarb/sqrt(n.lowcarb)), df = (n.lowcarb-1)), from = -20, to = 0, xname = 'mu')
for (mu in seq(from=-20, to= 0, by = 1)) {
  low.carb.t.test <- t.test(weightchange$weightchange[weightchange$diet=="Low Carb"], alternative = "greater", mu=mu)


abline(h=low.carb.t.test$p.value, v=mu)
}

Conf Curve

lowcarb.func <- function(mu) abs(2*(1-pt((xbar.lowcarb - mu)/(s.lowcarb/sqrt(n.lowcarb)), df = (n.lowcarb-1)))-1)
curve(lowcarb.func, from = -20, to = 0, xname = 'mu')

Checking Conf Curve

?t.test

low.carb.t.test <- t.test(weightchange$weightchange[weightchange$diet=="Low Carb"], conf.level = .80)

curve(lowcarb.func, from = -20, to = 0, xname = 'mu')
for (conf.int in seq(from=0, to= 1, by = .2)) {
  low.carb.t.test <- t.test(weightchange$weightchange[weightchange$diet=="Low Carb"], conf.level = conf.int)
  lowcarb.func <- function(mu) abs(2*(1-pt((xbar.lowcarb - mu)/(s.lowcarb/sqrt(n.lowcarb)), df = (n.lowcarb-1)))-1)

abline(h=conf.int, v=low.carb.t.test$conf.int)
}

Conf Density

lowcarb.dens.func <- function(mu) (1/(s.lowcarb/sqrt(n.lowcarb)))*(dt((xbar.lowcarb-mu)/((s.lowcarb/sqrt(n.lowcarb))), df=(n.lowcarb-1)))

curve(lowcarb.dens.func(mu), from = -20, to = 0, xname = "mu")

Checking Conf Density

?integrate
integrate(lowcarb.dens.func, lower = -Inf, upper = xbar.lowcarb)

curve(1-pt((xbar.lowcarb - mu)/(s.lowcarb/sqrt(n.lowcarb)), df = (n.lowcarb-1)), from = -20, to = 0, xname = 'mu')
for (mu in seq(from=-20, to= 0, by = .5)) {
  lowcarb.dens <- integrate(lowcarb.dens.func, lower = -Inf, upper = mu)$value


abline(h=lowcarb.dens, v=mu)
}

Low Fat

Conf Dist

xbar.lowfat <- mean(weightchange$weightchange[weightchange$diet=="Low Fat"])
s.lowfat <- sd(weightchange$weightchange[weightchange$diet=="Low Fat"])
n.lowfat <- length(weightchange$weightchange[weightchange$diet=="Low Fat"])
xbar.lowfat
s.lowfat
n.lowfat

curve(1-pt((xbar.lowfat - mu)/(s.lowfat/sqrt(n.lowfat)), df = (n.lowfat-1)), from = -20, to = 0, xname = 'mu')

Check

curve(1-pt((xbar.lowfat - mu)/(s.lowfat/sqrt(n.lowfat)), df = (n.lowfat-1)), from = -20, to = 0, xname = 'mu')
for (mu in seq(from=-20, to= 0, by = 1)) {
  low.fat.t.test <- t.test(weightchange$weightchange[weightchange$diet=="Low Fat"], alternative = "greater", mu=mu)


abline(h=low.fat.t.test$p.value, v=mu)
}

Conf Curve

lowfat.func <- function(mu) abs(2*(1-pt((xbar.lowfat - mu)/(s.lowfat/sqrt(n.lowfat)), df = (n.lowfat-1)))-1)
curve(lowfat.func, from = -20, to = 0, xname = 'mu')

low.fat.t.test <- t.test(weightchange$weightchange[weightchange$diet=="Low Fat"], conf.level = .80)

curve(lowfat.func, from = -20, to = 0, xname = 'mu')
for (conf.int in seq(from=0, to= 1, by = .2)) {
  low.fat.t.test <- t.test(weightchange$weightchange[weightchange$diet=="Low Fat"], conf.level = conf.int)
  lowcarb.func <- function(mu) abs(2*(1-pt((xbar.lowcarb - mu)/(s.lowcarb/sqrt(n.lowcarb)), df = (n.lowcarb-1)))-1)

abline(h=conf.int, v=low.fat.t.test$conf.int)
}

Conf Density

lowfat.dens.func <- function(mu) (1/(s.lowfat/sqrt(n.lowfat)))*(dt((xbar.lowfat-mu)/((s.lowfat/sqrt(n.lowfat))), df=(n.lowfat-1)))

curve(lowfat.dens.func(mu), from = -20, to = 0, xname = "mu")

Checking Conf Density

integrate(lowfat.dens.func, lower = -Inf, upper = xbar.lowfat)

curve(1-pt((xbar.lowfat - mu)/(s.lowfat/sqrt(n.lowfat)), df = (n.lowfat-1)), from = -20, to = 0, xname = 'mu')
for (mu in seq(from=-20, to= 0, by = .5)) {
  lowfat.dens <- integrate(lowfat.dens.func, lower = -Inf, upper = mu)$value


abline(h=lowfat.dens, v=mu)
}

Welch's t-test

weightchange.welch <- t.test(weightchange$weightchange ~ weightchange$diet)

Conf Dist

wc.df <- (((((s.lowcarb^2)/n.lowcarb)+(s.lowfat^2)/n.lowfat))^2)/((((s.lowcarb^2)/(n.lowcarb-1))^2)+(((s.lowfat^2)/(n.lowfat-1))^2))
wc.df # Weight Change Degrees of Freedom for Welch's t-test, should be 603.97 for this data

sample.diff <- xbar.lowcarb - xbar.lowfat

curve(1-pt((xbar.lowcarb - xbar.lowfat - delta)/(sqrt((((s.lowcarb)^2)/n.lowcarb)+((s.lowfat)^2)/n.lowfat)), df = 603.97), from = -10, to = 5, xname = 'delta')

Checking

curve(1-pt((xbar.lowcarb - xbar.lowfat - delta)/(sqrt((((s.lowcarb)^2)/n.lowcarb)+((s.lowfat)^2)/n.lowfat)), df = 603.97), from = -10, to = 5, xname = 'delta')
abline(h=.5, v=sample.diff)

Conf Curv

two.samp.func <- function(delta) abs(2*(1-pt((xbar.lowcarb - xbar.lowfat - delta)/(sqrt((((s.lowcarb)^2)/n.lowcarb)+((s.lowfat)^2)/n.lowfat)), df = 603.97))-1)

curve(two.samp.func, from = -10, to = 5, xname = 'delta')

Checking

curve(two.samp.func, from = -10, to = 5, xname = 'delta')
abline(h=.95, v = weightchange.welch$conf.int)

Conf Density

welch.dens.func <- function(delta) (1/sqrt((s.lowcarb^2/(n.lowcarb))+(s.lowfat^2/(n.lowfat)))*(dt((xbar.lowcarb- xbar.lowfat- delta)/sqrt((s.lowcarb^2/(n.lowcarb))+(s.lowfat^2/(n.lowfat))), df=(n.lowcarb-1))))

curve(lowcarb.dens.func(mu), from = -20, to = 0, xname = "mu")

Checking

integrate(welch.dens.func, lower = -Inf, upper = sample.diff)

curve(1-pt((xbar.lowcarb - xbar.lowfat - delta)/(sqrt((((s.lowcarb)^2)/n.lowcarb)+((s.lowfat)^2)/n.lowfat)), df = 603.97), from = -10, to = 5, xname = 'delta')
for (delta in seq(from=-20, to= 0, by = .5)) {
  welch.dens <- integrate(welch.dens.func, lower = -Inf, upper = delta)$value


abline(h=welch.dens, v=delta)
}

ddarmon/confdist documentation built on Aug. 8, 2020, 4:44 p.m.