# a.Create 10,000 simulations of a n = 100 sample.  Calculate the mean of the means and the standard deviation of the means

set.seed(123)

n1 = 100
n1_simulations = n2_simulations = 10000

# Simulate the means for n = 100
means_n1 = replicate(n1_simulations, mean(runif(n1, min = 3, max = 7)))

# Calculate the mean of the sample means and standard deviation
mean_of_means_n1 = mean(means_n1)
std_of_means_n1 = sd(means_n1)

cat("For n = 100:\n")

## For n = 100:

cat("Mean of the sample means:",mean_of_means_n1, "\n")

## Mean of the sample means: 4.998116

cat("Standard deviation of the sample means:", std_of_means_n1, "\n\n")

## Standard deviation of the sample means: 0.114493

# b.Create 10,000 simulations of a n = 10,000 sample. Calculate the mean of the means and the standard deviation of the means.

n2 = 10000
means_n2 = replicate(n2_simulations, mean(runif(n2, min = 3, max = 7)))

# Calculate the mean of the sample means and standard deviation
mean_of_means_n2 = mean(means_n2)
std_of_means_n2 = sd(means_n2)

# Print results for n = 10,000
cat("For n = 10,000:\n")

## For n = 10,000:

cat("Mean of the sample means:", mean_of_means_n2, "\n")

## Mean of the sample means: 4.999957

cat("Standard deviation of the sample means:", std_of_means_n2, "\n")

## Standard deviation of the sample means: 0.01155698

# 2.Now do the same with a distribution given as one normal distribution divided by another normal distribution.  What happens to the means and standard deviations?

set.seed(123)

# Function to generate one normal distribution divided by another normal distribution
generate_ratio_normal = function(n) {
  x = rnorm(n, mean = 5, sd = 1)  # Numerator normal distribution
  y = rnorm(n, mean = 2, sd = 1)  # Denominator normal distribution
  return(x / y)  # Ratio of the two normal distributions
}

n1 = 100
n1_simulations = n2_simulations = 10000

# Simulate the means for n = 100 using the ratio of two normal distributions
means_n1 = replicate(n1_simulations, mean(generate_ratio_normal(n1)))

# Calculate the mean of the sample means and standard deviation

mean_of_means_n1 = mean(means_n1)
sd_of_means_n1 = sd(means_n1)

cat("For n = 100 (Ratio of Normals):\n")

## For n = 100 (Ratio of Normals):

cat("Mean of the sample means:", mean_of_means_n1, "\n")

## Mean of the sample means: 3.690939

cat("Standard deviation of the sample means:", sd_of_means_n1, "\n\n")

## Standard deviation of the sample means: 36.65768

n2 = 10000

# Simulate the means for n = 10,000 using the ratio of two normal distributions
means_n2 = replicate(n2_simulations, mean(generate_ratio_normal(n2)))

# Calculate the mean of the sample means and standard deviation
mean_of_means_n2 = mean(means_n2)
std_of_means_n2 = sd(means_n2)

# Print results for n = 10,000
cat("For n = 10,000 (Ratio of Normals):\n")

## For n = 10,000 (Ratio of Normals):

cat("Mean of the sample means:", mean_of_means_n2, "\n")

## Mean of the sample means: 3.759293

cat("Standard deviation of the sample means:", std_of_means_n2, "\n")

## Standard deviation of the sample means: 75.03878