Statystyka/zajecia5/.Rhistory

load(url("http://ls.home.amu.edu.pl/data_sets/Centrala.RData"))
library(EnvStats)
epois(Centrala$Liczba,
method = "mle/mme/mvue",
ci = TRUE, ci.type = "two-sided", conf.level = 0.95,
ci.method = "exact")$interval$limits
epois(Centrala$Liczba,
method = "mle/mme/mvue",
ci = TRUE, ci.type = "two-sided", conf.level = 0.95,
ci.method = "pearson.hartley.approx")$interval$limits
epois(Centrala$Liczba,
method = "mle/mme/mvue",
ci = TRUE, ci.type = "two-sided", conf.level = 0.95,
ci.method = "normal.approx")$interval$limits
View(Centrala)
awarie <- read.table("http://ls.home.amu.edu.pl/data_sets/awarie.txt")
lambda_limits <- eexp(awarie$V1,
method = "mle/mme",
ci = TRUE, ci.type = "two-sided", conf.level = 0.95,
ci.method = "exact")$interval$limits
rev(1 / lambda_limits)
rev(1 / lambda_limits^2)
median_cint <- function(x, conf_level = 0.95) {
m <- mean(x^2)
a <- 1 - conf_level
n <- length(x)
L <- sqrt(log(2) * m * (1 - qnorm(1 - a / 2) / sqrt(n)))
R <- sqrt(log(2) * m * (1 + qnorm(1 - a / 2) / sqrt(n)))
w <- list(title = "mediana", est = sqrt(log(2) * m), l = L, r = R, conf_level = conf_level)
class(w) <- "confint"
return(w)
}
print.confint <- function(x) {
cat(x$conf_level * 100, "percent confidence interval:", "\n")
cat(x$l, " ", x$r, "\n")
}
summary.confint <- function(x) {
cat("\n", "Confidence interval of", x$title, "\n", "\n")
cat(x$conf_level * 100, "percent confidence interval:", "\n")
cat(x$l, " ", x$r, "\n")
cat("sample estimate", "\n")
cat(x$est, "\n")
}
x <- c(0.9, 6.2, 2.1, 4.1, 7.3,
1.0, 4.6, 6.4, 3.8, 5.0,
2.7, 9.2, 5.9, 7.4, 3.0,
4.9, 8.2, 5.0, 1.2, 10.1,
12.2, 2.8, 5.9, 8.2, 0.5)
median_cint(x)
summary(median_cint(x))
conf_int <- function(x, conf_level = 0.95) {
c(mean(x) - sd(x) * qt(1 - (1 - conf_level) / 2, length(x) - 1) / sqrt(length(x)),
mean(x) + sd(x) * qt(1 - (1 - conf_level) / 2, length(x) - 1) / sqrt(length(x)))
}
mu <- 1
sigma <- 3
nr <- 1000
df <- 3
lambda <- 3
n <- 10
cat(paste("n = ", n, sep = ""), "\n")
temp <- numeric(nr)
set.seed(1234)
for (i_nr in seq_len(nr)) {
x <- rnorm(n, mu, sigma)
temp_int <- conf_int(x)
temp[i_nr] <- ((temp_int[1] < mu) & (mu < temp_int[2]))
}
mean(temp)
temp <- numeric(nr)
set.seed(1234)
for (i_nr in seq_len(nr)) {
x <- rchisq(n, df)
temp_int <- conf_int(x)
temp[i_nr] <- ((temp_int[1] < df) & (df < temp_int[2]))
}
mean(temp)
temp <- numeric(nr)
set.seed(1234)
for (i_nr in seq_len(nr)) {
x <- rexp(n, rate = lambda)
temp_int <- conf_int(x)
temp[i_nr] <- ((temp_int[1] < 1 / lambda) & (1 / lambda < temp_int[2]))
}
mean(temp)
n <- 50
cat(paste("n = ", n, sep = ""), "\n")
temp <- numeric(nr)
set.seed(1234)
for (i_nr in seq_len(nr)) {
x <- rnorm(n, mu, sigma)
temp_int <- conf_int(x)
temp[i_nr] <- ((temp_int[1] < mu) & (mu < temp_int[2]))
}
mean(temp)
temp <- numeric(nr)
set.seed(1234)
for (i_nr in seq_len(nr)) {
x <- rchisq(n, df)
temp_int <- conf_int(x)
temp[i_nr] <- ((temp_int[1] < df) & (df < temp_int[2]))
}
mean(temp)
temp <- numeric(nr)
set.seed(1234)
for (i_nr in seq_len(nr)) {
x <- rexp(n, rate = lambda)
temp_int <- conf_int(x)
temp[i_nr] <- ((temp_int[1] < 1 / lambda) & (1 / lambda < temp_int[2]))
}
mean(temp)
n <- 100
cat(paste("n = ", n, sep = ""), "\n")
temp <- numeric(nr)
set.seed(1234)
for (i_nr in seq_len(nr)) {
x <- rnorm(n, mu, sigma)
temp_int <- conf_int(x)
temp[i_nr] <- ((temp_int[1] < mu) & (mu < temp_int[2]))
}
mean(temp)
temp <- numeric(nr)
set.seed(1234)
for (i_nr in seq_len(nr)) {
x <- rchisq(n, df)
temp_int <- conf_int(x)
temp[i_nr] <- ((temp_int[1] < df) & (df < temp_int[2]))
}
mean(temp)
temp <- numeric(nr)
set.seed(1234)
for (i_nr in seq_len(nr)) {
x <- rexp(n, rate = lambda)
temp_int <- conf_int(x)
temp[i_nr] <- ((temp_int[1] < 1 / lambda) & (1 / lambda < temp_int[2]))
}
mean(temp)