SimDesign in R

#install.packages("SimDesign")
library(SimDesign)
library(ggplot2)
set.seed(1)

# 1. DESIGN: Test different sample sizes (n)
Design <- createDesign(sample_size = c(10,100,1000, 10000, 100000))

# 2. GENERATE: Create random points
Generate <- function(condition, fixed_objects = NULL) {
  # Generate x and y coordinates between -1 and 1
  # condition$sample_size = the number of points for each trial
  x <- runif(condition$sample_size, min = -1, max = 1) 
  y <- runif(condition$sample_size, min = -1, max = 1)
  
  # Return the coordinates as a data frame
  data.frame(x = x, y = y)
}

# 3. ANALYSE: Calculate the estimate for pi
Analyse <- function(condition, dat, fixed_objects = NULL) {
  # radius of the circle is 1
  # points inside/on the circle has a distant to the origin <=1, so distance^2 <=1
  inside <- sum(dat$x^2 + dat$y^2 <= 1) # count how many points lie inside/on the circle
  
  # Pi approximation = 4 * (points_inside / total_points)
  pi_est <- 4 * (inside / condition$sample_size)
  
  # Return the single estimate
  return(pi_est)
}

# 4. SUMMARISE: Evaluate accuracy across replications
Summarise <- function(condition, results, fixed_objects = NULL) {
  # Calculate mean estimate and bias from true pi
  ret <- c(mean_pi = mean(results), 
           bias = mean(results) - pi,
           sd = sd(results))
  return(ret)
}

# 5. RUN SIMULATION
FinalResults <- runSimulation(design = Design, 
                             replications = 100, 
                             generate = Generate, 
                             analyse = Analyse, 
                             summarise = Summarise)

print(FinalResults)

set.seed(789)

library(tidyverse)
ames_housing <- read_csv("ames_housing.csv")

# We isolate the complete, clean rows where SalePrice and GrLivArea are not NA
ames_clean <- ames_housing[!is.na(ames_housing$SalePrice) & !is.na(ames_housing$GrLivArea), 
                       c("SalePrice", "GrLivArea")]


Design_Ames <- createDesign(impute_method = c("mean", "linear_regression", "hot_deck"))


Generate_Ames<- function(condition, fixed_objects) {
  
  df <- fixed_objects$clean_data
  n_rows <- nrow(df)
  
  # Keep a record of the true SalePrice values before we corrupt it
  true_SalePrice <- df$SalePrice
  
  # Intentionally punch missing holes (NA) into 20% of our real rows at random
  missing_indices <- sample(1:n_rows, size = round(0.20 * n_rows))
  df$SalePrice[missing_indices] <- NA
  
  return(list(corrupted_df = df, true_values = true_SalePrice, missing_ids = missing_indices))

}


Analyse_Ames <- function(condition, dat, fixed_objects) {
  test_df     <- dat$corrupted_df
  true_vector <- dat$true_values
  missing_ids <- dat$missing_ids
  
  # Track the true values of just the specific rows we deleted
  target_true <- true_vector[missing_ids]
  
  # Replace NA with the mean 
  if (condition$impute_method == "mean") {
    calculated_mean <- mean(test_df$SalePrice, na.rm = TRUE)
    imputed_values  <- rep(calculated_mean, length(missing_ids))
  }
  
  # Use a linear regression model to predict SalePrice from GrLivArea
  else if (condition$impute_method == "linear_regression") {
    fit            <- lm(SalePrice ~ GrLivArea, data = test_df)
    
    # Predict the missing SalePrice values using the corresponding GrLivArea values
    imputed_values <- predict(fit, newdata = test_df[missing_ids, ])
    
    # If regression outputs an NA due to missing GrLivArea, use overall mean
    imputed_values[is.na(imputed_values)] <- mean(test_df$SalePrice, na.rm = TRUE)
  }
  
  # Hot Deck Imputation: Randomly sample from rows that aren't missing and fill the missing values with them
  else if (condition$impute_method == "hot_deck") {
    valid_pool     <- test_df$SalePrice[!is.na(test_df$SalePrice)]
    imputed_values <- sample(valid_pool, size = length(missing_ids), replace = TRUE)
  }
  
  # Calculate estimation metrics for this run
  error <- imputed_values - target_true
  rmse  <- sqrt(mean(error^2))
  mean_imputed <- mean(imputed_values)
  mean_true    <- mean(target_true)
  
  return(c(RMSE = rmse, Imputed_Mean = mean_imputed, True_Mean = mean_true))
}


Summarise_Ames <- function(condition, results, fixed_objects) {
  bias_score <- mean(results[, "Imputed_Mean"]) - mean(results[, "True_Mean"])
  
  ret <- c(
    Average_RMSE = mean(results[, "RMSE"]),
    Empirical_Bias = bias_score
  )
  return(ret)
}


SimulationResults_Ames <- runSimulation(design = Design_Ames,
                                       replications = 200,
                                       generate = Generate_Ames,
                                       analyse = Analyse_Ames,
                                       summarise = Summarise_Ames,
                                       fixed_objects = list(clean_data = ames_clean))

print(SimulationResults_Ames)