Solved problem 9.

This commit is contained in:
Filipe Rodrigues 2023-06-16 14:04:56 +01:00
parent 12627af349
commit 149b46c5a1
5 changed files with 80 additions and 0 deletions

BIN
9/P9_1.png Normal file

Binary file not shown.

After

Width:  |  Height:  |  Size: 54 KiB

BIN
9/P9_2.png Normal file

Binary file not shown.

After

Width:  |  Height:  |  Size: 52 KiB

BIN
9/P9_3.png Normal file

Binary file not shown.

After

Width:  |  Height:  |  Size: 47 KiB

76
9/script.R Executable file
View File

@ -0,0 +1,76 @@
#!/bin/env Rscript
# Libraries
library("pracma")
library("ggplot2")
# Set the seed for the RNG
set.seed(1549)
# All the `n`s we'll be iterating over
all_n <- list(30, 50, 100, 200, 300, 500, 1000)
all_n <- sapply(all_n, function(n) n)
print(all_n)
# Calculate the samples
p <- 0.8
samples_len <- 2000
samples <- lapply(all_n, function(n) {
lapply(1:samples_len, function(sample_idx) {
rbinom(n, size = 1, prob = p)
})
})
# Calculate the interval lengths for each method
confidence <- 0.98
z <- qnorm((1.0 + confidence) / 2.0)
method1_length <- function(mean, n) {
a <- z^2 / n
b <- a * (a + 4 * mean - 4 * mean^2)
p_min <- (2 * mean + a - sqrt(b)) / (2 * (a + 1))
p_max <- (2 * mean + a + sqrt(b)) / (2 * (a + 1))
p_max - p_min
}
method2_length <- function(mean, n) {
a <- z^2 * (mean * (1 - mean)) / n
p_min <- -sqrt(a) + mean
p_max <- sqrt(a) + mean
p_max - p_min
}
lengths1 <- mapply(function(samples_n, n) {
mean(sapply(samples_n, function(sample) {
method1_length(mean(sample), n)
}))
}, samples, all_n)
print(lengths1)
lengths2 <- mapply(function(samples_n, n) {
mean(sapply(samples_n, function(sample) {
method2_length(mean(sample), n)
}))
}, samples, all_n)
print(lengths2)
# Finally plot and save it
df <- data.frame(
all_n = all_n,
lengths1 = lengths1,
lengths2 = lengths2
)
plot <- ggplot(df, aes(x = all_n)) +
geom_line(aes(y = lengths1, color = "Method 1")) +
geom_line(aes(y = lengths2, color = "Method 2")) +
xlab("N") +
ylab("") +
scale_color_manual(values = c("red", "blue")) +
theme(legend.title = element_blank())
ggsave(plot = plot, "output.png")

View File

@ -16,6 +16,10 @@ https://docs.google.com/document/d/1DnvOxo3GlUY7sBHUCcyMEbMlQE2AYFBWAOSrNZyb4d4
https://docs.google.com/document/d/1x_AsfHc4JNxV8P6oOfUoZCu1yLQ9TQ82zRSoT7pDRrM
## 9
https://docs.google.com/document/d/1TmL1dgG8DoG6-_J2rQPIU_SOcNpKvbDBL0KVmuLRe1g
# Util
## 6