Direct-Mail Fundraising Optimization

Veterans’ Organization Case Study

Olivia Shipley

2025-08-02

knitr::opts_chunk$set(echo = TRUE, warning = FALSE, message = FALSE)
library(tidyverse)
library(randomForest)
library(tidyr)
library(caret)
library(corrplot)
library(VIM)
library(ROCR)
library(pROC)
library(gbm)
library(conflicted)
library(dplyr)

conflicts_prefer(dplyr::filter)
conflicts_prefer(dplyr::select)

set.seed(12345)

A Data-Driven Transformation Story

This analysis follows the journey of a veterans’ organization facing a classic business challenge: how to raise funds efficiently without wasting resources. Using Disney’s Story Spine structure, our narrative unfolds as:

This analysis transforms a veterans’ organization’s fundraising approach from costly mass mailings to precision-targeted campaigns using predictive modeling. The current mass mailing strategy loses approximately $17 per 1,000 mailings, while our data-driven approach projects positive ROI through strategic donor identification.

🎯 Business Objectives & Goals

Primary Goal

Transform fundraising operations from cost-ineffective mass mailings to profitable, targeted donor engagement that maximizes resources available for veteran services.

Specific Business Objectives

  1. Cost Reduction: Eliminate wasteful mailings to unlikely donors (currently 94.9% of mass mailings)
  2. Revenue Optimization: Focus resources on high-probability donors to maximize return on investment
  3. Scalability: Develop a repeatable, data-driven process applicable to the organization’s 13-million donor database
  4. Mission Alignment: Increase funds available for veteran services through improved operational efficiency

Success Metrics

📊 Business Problem Analysis

#current funraising analysis
response_rate <- 0.051  # overall response rate
avg_donation <- 13.00   # avg donation amount
cost_per_mail <- 0.68   # cost per mailing

#per 1000 mailings
mailings_per_batch <- 1000
expected_donors <- mailings_per_batch * response_rate
current_revenue <- expected_donors * avg_donation
current_cost <- mailings_per_batch * cost_per_mail
current_net <- current_revenue - current_cost
cat(sprintf("Current Mass Mailing Economics (per 1,000 mailings):\n"))
## Current Mass Mailing Economics (per 1,000 mailings):
cat(sprintf("  Expected donors: %.0f (%.1f%% response rate)\n", expected_donors, response_rate * 100))
##   Expected donors: 51 (5.1% response rate)
cat(sprintf("  Revenue: $%.2f\n", current_revenue))
##   Revenue: $663.00
cat(sprintf("  Total cost: $%.2f\n", current_cost))
##   Total cost: $680.00
cat(sprintf("  Net result: $%.2f\n", current_net))
##   Net result: $-17.00
cat(sprintf("  Status: %s\n\n", ifelse(current_net > 0, "Profitable", "LOSING MONEY")))
##   Status: LOSING MONEY

Objective: Develop a predictive model to identify likely donors

Goal: Maximize net profit by reducing mailings to non-donors

Success Metric: Net profit improvement over mass mailing approach

🗄️ Data Sources & Data Description

Data Acquisition

I was provided two datasets to complete this analysis:

  1. Training Dataset: 3,000 records with donor outcomes (STA6543 Summer2025 Fundraising.csv)
  2. Test Dataset: 120 records for prediction (Future Fundraising Summer 2025.csv)

Sampling Methodology: Weighted Sampling Justification

Critical Decision: Intially, a 50/50 weighted sampling was used, instead of simple random sampling.

Rationale:

Statistical Consideration: This sampling approach requires careful interpretation during model deployment, as we must account for the artificial balance when applying models to naturally imbalanced populations.

Variable Descriptions & Business Inputs

The dataset contains 20 predictor variables across four categories:

Geographic Variables: - zipconvert2-5: Zip code groupings (00000-19999, 20000-39999, 40000-59999, 60000-79999, 80000-99999)

Demographic Variables:

- homeowner: Property ownership status

- num_child: Number of children

- female: Gender indicator

- income: Household income level (1-7 scale)

Wealth Indicators:

- wealth: Relative wealth rating (0-9 scale)

- home_value: Average neighborhood home value (hundreds of dollars)

- med_fam_inc: Median family income in area (hundreds of dollars)

- avg_fam_inc: Average family income in area (hundreds of dollars)

- pct_lt15k: Percentage earning less than $15K in area

Behavioral/Engagement Variables:

- num_prom: Lifetime number of promotions received

- lifetime_gifts: Total dollar amount of historical donations

- largest_gift: Largest single donation amount

- last_gift: Most recent donation amount

- months_since_donate: Time since last donation (recency)

- time_lag: Time between first and second gifts

- avg_gift: Average donation amount per gift

🔬 Exploratory Data Analysis

The Challenge: Beyond Simple Demographics

Our initial hypothesis was that basic demographics (geography, homeowner status, gender) might drive donation patterns. However, as we’ll see, the story is more complex and requires more sophisticated behavioral modeling.

Methodology Selection and Rationale

Primary Approach: Binary classification with business-focused evaluation metrics

Why This Methodology:

  1. Business Alignment: The problem is fundamentally binary (mail/don’t mail)
  2. Cost-Benefit Focus: Evaluation prioritizes profit over accuracy
  3. Interpretability: Models must be explainable to stakeholders
  4. Scalability: Solution must work across large donor databases

Alternative Methodologies Considered:

train_data <- read.csv("STA6543 Summer2025 Fundraising.csv", stringsAsFactors = TRUE)
test_data <- read.csv("Future Fundraising Summer 2025.csv", stringsAsFactors = TRUE)
cat("Training Data Dimensions:", dim(train_data), "\n")
## Training Data Dimensions: 3000 21
cat("Test Data Dimensions:", dim(test_data), "\n")
## Test Data Dimensions: 120 20
table(train_data$target)
## 
##    Donor No Donor 
##     1499     1501
prop.table(table(train_data$target))
## 
##     Donor  No Donor 
## 0.4996667 0.5003333
str(train_data)
## 'data.frame':    3000 obs. of  21 variables:
##  $ zipconvert2        : Factor w/ 2 levels "No","Yes": 2 1 1 1 1 1 1 2 1 2 ...
##  $ zipconvert3        : Factor w/ 2 levels "No","Yes": 1 1 1 2 2 1 1 1 1 1 ...
##  $ zipconvert4        : Factor w/ 2 levels "No","Yes": 1 1 1 1 1 1 2 1 1 1 ...
##  $ zipconvert5        : Factor w/ 2 levels "No","Yes": 1 2 2 1 1 2 1 1 2 1 ...
##  $ homeowner          : Factor w/ 2 levels "No","Yes": 2 1 2 2 2 2 2 2 2 2 ...
##  $ num_child          : int  1 2 1 1 1 1 1 1 1 1 ...
##  $ income             : int  1 5 3 4 4 4 4 4 4 1 ...
##  $ female             : Factor w/ 2 levels "No","Yes": 1 2 1 1 2 2 1 2 2 2 ...
##  $ wealth             : int  7 8 4 8 8 8 5 8 8 5 ...
##  $ home_value         : int  698 828 1471 547 482 857 505 1438 1316 428 ...
##  $ med_fam_inc        : int  422 358 484 386 242 450 333 458 541 203 ...
##  $ avg_fam_inc        : int  463 376 546 432 275 498 388 533 575 271 ...
##  $ pct_lt15k          : int  4 13 4 7 28 5 16 8 11 39 ...
##  $ num_prom           : int  46 32 94 20 38 47 51 21 66 73 ...
##  $ lifetime_gifts     : num  94 30 177 23 73 139 63 26 108 161 ...
##  $ largest_gift       : num  12 10 10 11 10 20 15 16 12 6 ...
##  $ last_gift          : num  12 5 8 11 10 20 10 16 7 3 ...
##  $ months_since_donate: int  34 29 30 30 31 37 37 30 31 32 ...
##  $ time_lag           : int  6 7 3 6 3 3 8 6 1 7 ...
##  $ avg_gift           : num  9.4 4.29 7.08 7.67 7.3 ...
##  $ target             : Factor w/ 2 levels "Donor","No Donor": 1 1 2 2 1 1 1 2 1 1 ...
summary(train_data)
##  zipconvert2 zipconvert3 zipconvert4 zipconvert5 homeowner    num_child    
##  No :2352    No :2449    No :2357    No :1846    No : 688   Min.   :1.000  
##  Yes: 648    Yes: 551    Yes: 643    Yes:1154    Yes:2312   1st Qu.:1.000  
##                                                             Median :1.000  
##                                                             Mean   :1.069  
##                                                             3rd Qu.:1.000  
##                                                             Max.   :5.000  
##      income      female         wealth        home_value      med_fam_inc    
##  Min.   :1.000   No :1169   Min.   :0.000   Min.   :   0.0   Min.   :   0.0  
##  1st Qu.:3.000   Yes:1831   1st Qu.:5.000   1st Qu.: 554.8   1st Qu.: 278.0  
##  Median :4.000              Median :8.000   Median : 816.5   Median : 355.0  
##  Mean   :3.899              Mean   :6.396   Mean   :1143.3   Mean   : 388.4  
##  3rd Qu.:5.000              3rd Qu.:8.000   3rd Qu.:1341.2   3rd Qu.: 465.0  
##  Max.   :7.000              Max.   :9.000   Max.   :5945.0   Max.   :1500.0  
##   avg_fam_inc       pct_lt15k        num_prom      lifetime_gifts  
##  Min.   :   0.0   Min.   : 0.00   Min.   : 11.00   Min.   :  15.0  
##  1st Qu.: 318.0   1st Qu.: 5.00   1st Qu.: 29.00   1st Qu.:  45.0  
##  Median : 396.0   Median :12.00   Median : 48.00   Median :  81.0  
##  Mean   : 432.3   Mean   :14.71   Mean   : 49.14   Mean   : 110.7  
##  3rd Qu.: 516.0   3rd Qu.:21.00   3rd Qu.: 65.00   3rd Qu.: 135.0  
##  Max.   :1331.0   Max.   :90.00   Max.   :157.00   Max.   :5674.9  
##   largest_gift       last_gift      months_since_donate    time_lag     
##  Min.   :   5.00   Min.   :  0.00   Min.   :17.00       Min.   : 0.000  
##  1st Qu.:  10.00   1st Qu.:  7.00   1st Qu.:29.00       1st Qu.: 3.000  
##  Median :  15.00   Median : 10.00   Median :31.00       Median : 5.000  
##  Mean   :  16.65   Mean   : 13.48   Mean   :31.13       Mean   : 6.876  
##  3rd Qu.:  20.00   3rd Qu.: 16.00   3rd Qu.:34.00       3rd Qu.: 9.000  
##  Max.   :1000.00   Max.   :219.00   Max.   :37.00       Max.   :77.000  
##     avg_gift            target    
##  Min.   :  2.139   Donor   :1499  
##  1st Qu.:  6.333   No Donor:1501  
##  Median :  9.000                  
##  Mean   : 10.669                  
##  3rd Qu.: 12.800                  
##  Max.   :122.167

Pattern Analysis

#num
numeric_vars <- train_data |> 
  select_if(is.numeric) |>
  select(-matches("zipconvert")) 
cor_matrix <- cor(numeric_vars, use = "complete.obs")

corrplot(cor_matrix, method = "color", type = "upper", 
         order = "hclust", addCoef.col = "black", 
         tl.cex = 0.8, number.cex = 0.7)

#correlated pairs
high_cor <- which(abs(cor_matrix) > 0.7 & abs(cor_matrix) < 1, arr.ind = TRUE)

if(nrow(high_cor) > 0) {
  cat("Highly correlated variable pairs (|r| > 0.7):\n")
  for(i in 1:nrow(high_cor)) {
    var1 <- rownames(cor_matrix)[high_cor[i,1]]
    var2 <- colnames(cor_matrix)[high_cor[i,2]]
    correlation <- cor_matrix[high_cor[i,1], high_cor[i,2]]
    cat(sprintf("  %s - %s: %.3f\n", var1, var2, correlation))
  }
} else {
  cat("No highly correlated variable pairs found (|r| > 0.7)\n\n")
}
## Highly correlated variable pairs (|r| > 0.7):
##   med_fam_inc - home_value: 0.738
##   avg_fam_inc - home_value: 0.753
##   home_value - med_fam_inc: 0.738
##   avg_fam_inc - med_fam_inc: 0.972
##   home_value - avg_fam_inc: 0.753
##   med_fam_inc - avg_fam_inc: 0.972
##   avg_gift - last_gift: 0.866
##   last_gift - avg_gift: 0.866
business_critical_vars <- c("lifetime_gifts", "largest_gift", "avg_gift","num_prom", "months_since_donate", "wealth", "home_value", "income")

business_critical_vars <- business_critical_vars[business_critical_vars %in% names(train_data)]

cat("Analyzing donor vs non-donor differences:\n\n")
## Analyzing donor vs non-donor differences:
for(var in business_critical_vars) {
  donor_avg <- mean(train_data[train_data$target == "Donor", var], na.rm = TRUE)
  non_donor_avg <- mean(train_data[train_data$target == "No Donor", var], na.rm = TRUE)
  difference <- donor_avg - non_donor_avg
  percent_diff <- (difference / non_donor_avg) * 100
  
  cat(sprintf("%s:\n", var))
  cat(sprintf("   Donors: %.2f | Non-donors: %.2f | Difference: %.1f%%\n", 
              donor_avg, non_donor_avg, percent_diff))

  #cute sig test
  t_test <- t.test(train_data[train_data$target == "Donor", var],
                   train_data[train_data$target == "No Donor", var])
  
  if(t_test$p.value < 0.001) {
    cat("   ⭐ HIGHLY SIGNIFICANT for prediction\n")
  } else if(t_test$p.value < 0.05) {
    cat("   ✅ Significant for prediction\n")
  } else {
    cat("   ❌ Not statistically reliable\n")
  }
  cat("\n")
  }
## lifetime_gifts:
##    Donors: 113.67 | Non-donors: 107.81 | Difference: 5.4%
##    ❌ Not statistically reliable
## 
## largest_gift:
##    Donors: 16.25 | Non-donors: 17.05 | Difference: -4.7%
##    ❌ Not statistically reliable
## 
## avg_gift:
##    Donors: 10.10 | Non-donors: 11.23 | Difference: -10.0%
##    ⭐ HIGHLY SIGNIFICANT for prediction
## 
## num_prom:
##    Donors: 50.70 | Non-donors: 47.58 | Difference: 6.5%
##    ⭐ HIGHLY SIGNIFICANT for prediction
## 
## months_since_donate:
##    Donors: 30.58 | Non-donors: 31.68 | Difference: -3.5%
##    ⭐ HIGHLY SIGNIFICANT for prediction
## 
## wealth:
##    Donors: 6.40 | Non-donors: 6.39 | Difference: 0.2%
##    ❌ Not statistically reliable
## 
## home_value:
##    Donors: 1163.82 | Non-donors: 1122.75 | Difference: 3.7%
##    ❌ Not statistically reliable
## 
## income:
##    Donors: 3.96 | Non-donors: 3.84 | Difference: 3.1%
##    ✅ Significant for prediction

The correlation analysis reveals significant multicollinearity between income-related variables (med_fam_inc and avg_fam_inc with r = 0.972) and between giving variables (avg_gift and last_gift with r = 0.866), which may require attention in linear models.

Surprisingly, the strongest predictors are behavioral rather than financial: num_prom (engagement), months_since_donate (recency), and avg_gift all show highly significant differences between donors and non-donors. Interestingly, donors actually have lower average gift amounts, suggesting frequent small donors rather than occasional large donors drive the organization’s fundraising success.

Donor Analysis

#zip codes
zip_summary <- train_data |>
  select(zipconvert2:zipconvert5, target) |>
  gather(zip_group, value, zipconvert2:zipconvert5) |>
  filter(value == "Yes") |>
  count(zip_group, target) |>
  group_by(zip_group) |>
  mutate(prop_donor = n/sum(n)) |>
  filter(target == "Donor")

cat("Donor rates by zip code group:\n")
## Donor rates by zip code group:
print(zip_summary)
## # A tibble: 4 × 4
## # Groups:   zip_group [4]
##   zip_group   target     n prop_donor
##   <chr>       <fct>  <int>      <dbl>
## 1 zipconvert2 Donor    320      0.494
## 2 zipconvert3 Donor    269      0.488
## 3 zipconvert4 Donor    318      0.495
## 4 zipconvert5 Donor    592      0.513

The highest donor rate, according to zip code is the “zipconver5” group with 51.3%. The lowest is “zipconvert3” at 49.8%. With only a 2.5% difference between the highest and lowest zip areas, it does not seem like geography is a major predictor in this case.

#demo (homeowner & gender)
cat_vars <- c("homeowner", "female")

for(var in cat_vars) {
  cat(sprintf("\n%s Analysis:", toupper(var)))
  tab <- table(train_data[[var]], train_data$target)
  print(tab)
  cat("Proportions:")
  print(round(prop.table(tab, 1), 3))

  #sig test  
  chi_test <- chisq.test(tab)
  cat(sprintf("Chi-square p-value: %.4f", chi_test$p.value))
  if(chi_test$p.value < 0.05) cat(" (Statistically significant)")
  cat("\n")
}
## 
## HOMEOWNER Analysis:     
##       Donor No Donor
##   No    327      361
##   Yes  1172     1140
## Proportions:     
##       Donor No Donor
##   No  0.475    0.525
##   Yes 0.507    0.493
## Chi-square p-value: 0.1576
## 
## FEMALE Analysis:     
##       Donor No Donor
##   No    566      603
##   Yes   933      898
## Proportions:     
##       Donor No Donor
##   No  0.484    0.516
##   Yes 0.510    0.490
## Chi-square p-value: 0.1873

The analysis of homeowner status show a chi-square p-value of 0.1576. Indicating a slightly positive trend but is not statistically significant. In terms of gender, self-identifying females show slightly higher giving rates (51.0%) compared to self-identifying males (48.4%). Again, with a chi-square p-value of 0.1873, the difference is not statistically significant.

#donor vs non-donor
key_business_vars <- c("lifetime_gifts", "largest_gift", "avg_gift", 
                      "num_prom", "months_since_donate", "wealth")

for(var in key_business_vars) {
  donor_data <- train_data[train_data$target == "Donor", var]
  non_donor_data <- train_data[train_data$target == "No Donor", var]
  
  t_test_result <- t.test(donor_data, non_donor_data)
  
  #cohen's d
  pooled_sd <- sqrt(((length(donor_data)-1) * var(donor_data) + 
                     (length(non_donor_data)-1) * var(non_donor_data)) / 
                    (length(donor_data) + length(non_donor_data) - 2))
  cohens_d <- (mean(donor_data) - mean(non_donor_data)) / pooled_sd
  
  cat(sprintf("\n%s:\n", var))
  cat(sprintf("  Donor average: %.2f\n", mean(donor_data)))
  cat(sprintf("  Non-donor average: %.2f\n", mean(non_donor_data)))
  cat(sprintf("  Difference: %.2f\n", mean(donor_data) - mean(non_donor_data)))
  cat(sprintf("  p-value: %.2e\n", t_test_result$p.value))
  cat(sprintf("  Effect size (Cohen's d): %.3f", cohens_d))
  
  # Interpret effect size
  if(abs(cohens_d) < 0.2) effect <- "(negligible)"
  else if(abs(cohens_d) < 0.5) effect <- "(small effect)"
  else if(abs(cohens_d) < 0.8) effect <- "(medium effect)"
  else effect <- "(large effect)"
  cat(sprintf(" %s\n", effect))
}
## 
## lifetime_gifts:
##   Donor average: 113.67
##   Non-donor average: 107.81
##   Difference: 5.86
##   p-value: 2.82e-01
##   Effect size (Cohen's d): 0.039 (negligible)
## 
## largest_gift:
##   Donor average: 16.25
##   Non-donor average: 17.05
##   Difference: -0.80
##   p-value: 3.30e-01
##   Effect size (Cohen's d): -0.036 (negligible)
## 
## avg_gift:
##   Donor average: 10.10
##   Non-donor average: 11.23
##   Difference: -1.13
##   p-value: 3.34e-05
##   Effect size (Cohen's d): -0.152 (negligible)
## 
## num_prom:
##   Donor average: 50.70
##   Non-donor average: 47.58
##   Difference: 3.11
##   p-value: 1.79e-04
##   Effect size (Cohen's d): 0.137 (negligible)
## 
## months_since_donate:
##   Donor average: 30.58
##   Non-donor average: 31.68
##   Difference: -1.10
##   p-value: 1.86e-13
##   Effect size (Cohen's d): -0.270 (small effect)
## 
## wealth:
##   Donor average: 6.40
##   Non-donor average: 6.39
##   Difference: 0.02
##   p-value: 8.65e-01
##   Effect size (Cohen's d): 0.006 (negligible)

⚙️ Variable Transformations & Feature Engineering

Exclusions Applied

No variable exclusions were implemented for the following reasons:

1. Domain Knowledge: All variables represent potentially valuable donor insights 2. Model Robustness: Random Forest can handle correlated variables effectively 3. Business Requirements: Stakeholders requested comprehensive analysis of all available data

Future Consideration: In production deployment, variables with consistently low importance scores across multiple model iterations could be excluded to simplify the model.

Variable Transformations

Categorical Encoding: All categorical variables were converted to factors for proper handling in R models.

No Additional Transformations Applied:

📈 Model Development

Data Partitioning

#80/20 train/validation split
train_index <- createDataPartition(train_data$target, p = 0.8, list = FALSE)
train_set <- train_data[train_index, ]
validation_set <- train_data[-train_index, ]

cat("*Data Partitioning (80/20 split):\n")
## *Data Partitioning (80/20 split):
cat(sprintf("  Training set: %d observations\n", nrow(train_set)))
##   Training set: 2401 observations
cat(sprintf("  Validation set: %d observations\n", nrow(validation_set)))
##   Validation set: 599 observations
#verify target dist
cat("\n*Target distribution verification:\n")
## 
## *Target distribution verification:
cat("Training set:")
## Training set:
print(table(train_set$target))
## 
##    Donor No Donor 
##     1200     1201
cat("Validation set:")
## Validation set:
print(table(validation_set$target))
## 
##    Donor No Donor 
##      299      300

Model Selection Strategy

I chose to build two complementary models:

Logistic Regression: Provides interpretable coefficients that executives can understand and trust. The simplicity makes it easy to explain why certain donors are targeted.

Random Forest: Captures complex interactions between variables that logistic regression might miss. Particularly valuable for understanding non-linear relationships in donor behavior.

Both models will be evaluated primarily on net profit rather than accuracy, since our business objective is to maximize fundraising efficiency, not just prediction accuracy.

Model 1: Logistic Regression

Rationale: Provides interpretable coefficients for executive presentations and regulatory compliance.

Benefits:

- Transparency in decision-making

- Fast training and prediction

- Probability estimates for threshold optimization

- Well-understood statistical properties

log_model <- glm(target ~ ., data = train_set, family = "binomial")

cat("\nModel Summary (Top 10 Most Significant Variables):\n")
## 
## Model Summary (Top 10 Most Significant Variables):
coef_summary <- summary(log_model)$coefficients
significant_vars <- coef_summary[order(coef_summary[,4])[1:10], ]
print(round(significant_vars, 4))
##                     Estimate Std. Error z value Pr(>|z|)
## months_since_donate   0.0689     0.0114  6.0531   0.0000
## income               -0.0831     0.0293 -2.8378   0.0045
## num_child             0.2681     0.1262  2.1238   0.0337
## home_value           -0.0002     0.0001 -2.0803   0.0375
## avg_fam_inc           0.0017     0.0011  1.5023   0.1330
## num_prom             -0.0035     0.0027 -1.3126   0.1893
## last_gift             0.0097     0.0090  1.0844   0.2782
## med_fam_inc          -0.0011     0.0010 -1.0402   0.2982
## wealth               -0.0167     0.0203 -0.8242   0.4098
## homeownerYes         -0.0688     0.1056 -0.6510   0.5151
log_pred_prob <- predict(log_model, validation_set, type = "response")
log_pred_class <- ifelse(log_pred_prob > 0.5, "Donor", "No Donor")

# Calculate business impact for Logistic Regression
log_tp <- sum(log_pred_class == "Donor" & validation_set$target == "Donor")
log_fp <- sum(log_pred_class == "Donor" & validation_set$target == "No Donor")
log_fn <- sum(log_pred_class == "No Donor" & validation_set$target == "Donor")
log_profit <- log_tp * 13.00 - (log_tp + log_fp) * 0.68
log_accuracy <- mean(log_pred_class == validation_set$target)

cat("Logistic Regression Results:\n")
## Logistic Regression Results:
cat(sprintf("Accuracy: %.1f%%, Net Profit: $%.2f\n", log_accuracy * 100, log_profit))
## Accuracy: 45.4%, Net Profit: $1435.96

Model 2: Random Forest

Rationale: Captures non-linear relationships and variable interactions without manual feature engineering.

Benefits: - Robust to outliers and missing values

- Handles mixed data types effectively

- Provides variable importance rankings

- Excellent performance on tabular data

rf_model <- randomForest(target ~ ., data = train_set, 
                        ntree = 500, 
                        mtry = sqrt(ncol(train_set) - 1),
                        importance = TRUE)

cat(sprintf("\nModel Configuration:\n"))
## 
## Model Configuration:
cat(sprintf("  Number of trees: %d\n", rf_model$ntree))
##   Number of trees: 500
cat(sprintf("  Variables per split: %d\n", rf_model$mtry))
##   Variables per split: 4
cat(sprintf("  OOB Error Rate: %.2f%%\n", rf_model$err.rate[nrow(rf_model$err.rate), "OOB"] * 100))
##   OOB Error Rate: 46.56%
cat("\nTop 10 Most Important Variables:\n")
## 
## Top 10 Most Important Variables:
importance_scores <- importance(rf_model)
var_importance <- data.frame(
  Variable = rownames(importance_scores),
  MeanDecreaseGini = importance_scores[, "MeanDecreaseGini"]
) |>
  arrange(desc(MeanDecreaseGini)) |>
  head(10)
print(var_importance)
##                                Variable MeanDecreaseGini
## home_value                   home_value        116.91196
## avg_gift                       avg_gift        112.61514
## med_fam_inc                 med_fam_inc        107.48970
## avg_fam_inc                 avg_fam_inc        105.64319
## lifetime_gifts           lifetime_gifts        105.20543
## num_prom                       num_prom         97.23649
## pct_lt15k                     pct_lt15k         90.42228
## time_lag                       time_lag         78.20972
## months_since_donate months_since_donate         76.11002
## last_gift                     last_gift         65.04783
rf_pred_class <- predict(rf_model, validation_set)
rf_pred_prob <- predict(rf_model, validation_set, type = "prob")[, "Donor"]

# Calculate business impact for Random Forest
rf_tp <- sum(rf_pred_class == "Donor" & validation_set$target == "Donor")
rf_fp <- sum(rf_pred_class == "Donor" & validation_set$target == "No Donor")
rf_fn <- sum(rf_pred_class == "No Donor" & validation_set$target == "Donor")
rf_profit <- rf_tp * 13.00 - (rf_tp + rf_fp) * 0.68
rf_accuracy <- mean(rf_pred_class == validation_set$target)

cat("Random Forest Results:\n")
## Random Forest Results:
cat(sprintf("Accuracy: %.1f%%, Net Profit: $%.2f\n", rf_accuracy * 100, rf_profit))
## Accuracy: 51.8%, Net Profit: $1915.76

To ensure robust model performance, we’ll conduct cross-validation analysis:

⚠️⚠️⚠️ Cross-Validation Analysis

train_set$target <- factor(train_set$target, 
                          levels = c("No Donor", "Donor"),
                          labels = c("NoDonor", "Donor"))

validation_set$target <- factor(validation_set$target,
                               levels = c("No Donor", "Donor"), 
                               labels = c("NoDonor", "Donor"))
library(caret)

ctrl <- trainControl(
  method = "cv", 
  number = 5,
  summaryFunction = twoClassSummary,
  classProbs = TRUE,
  savePredictions = TRUE
)

#cv for rf
set.seed(12345)
rf_cv <- train(target ~ ., 
               data = train_set,
               method = "rf",
               trControl = ctrl,
               metric = "ROC",
               ntree = 500)

#cv for logre
log_cv <- train(target ~ ., 
                data = train_set,
                method = "glm",
                family = "binomial",
                trControl = ctrl,
                metric = "ROC")

cat("Cross-Validation Results:\n")
## Cross-Validation Results:
cat(sprintf("Random Forest CV ROC: %.3f (±%.3f)\n", 
            mean(rf_cv$resample$ROC), sd(rf_cv$resample$ROC)))
## Random Forest CV ROC: 0.569 (±0.020)
cat(sprintf("Logistic Regression CV ROC: %.3f (±%.3f)\n", 
            mean(log_cv$resample$ROC), sd(log_cv$resample$ROC)))
## Logistic Regression CV ROC: 0.584 (±0.019)
cv_business_results <- rf_cv$pred |>
  group_by(Resample) |>
  summarise(
    business_profit = sum(ifelse(pred == "Donor" & obs == "Donor", 13.00 - 0.68,
                                ifelse(pred == "Donor" & obs == "No Donor", -0.68, 0)))
  )

cat(sprintf("Random Forest CV Business Profit: $%.2f (±%.2f)\n", 
            mean(cv_business_results$business_profit), 
            sd(cv_business_results$business_profit)))
## Random Forest CV Business Profit: $4740.74 (±288.98)
#convert validation predictions back
validation_set_original <- validation_set
validation_set_original$target <- factor(validation_set$target,
                                         levels = c("NoDonor", "Donor"),
                                         labels = c("No Donor", "Donor"))

log_pred_class_original <- factor(log_pred_class, 
                                 levels = c("Donor", "No Donor"))
rf_pred_class_original <- factor(rf_pred_class,
                                levels = c("Donor", "No Donor"))

# Recalculate for comparison table
log_performance <- list(net_profit =log_profit, accuracy =log_accuracy,                      precision = log_tp/(log_tp + log_fp), 
                     recall = log_tp/(log_tp + log_fn))

rf_performance <- list(net_profit = rf_profit, accuracy = rf_accuracy,
                      precision = rf_tp/(rf_tp + rf_fp), 
                      recall = rf_tp/(rf_tp + rf_fn))


model_comparison <- data.frame(
  Model = c("Logistic Regression", "Random Forest"),
  Net_Profit = c(log_profit, rf_profit),
  Accuracy = c(log_accuracy, rf_accuracy) * 100,
  Precision = c(log_tp/(log_tp + log_fp), rf_tp/(rf_tp + rf_fp)) * 100,
  Recall = c(log_tp/(log_tp + log_fn), rf_tp/(rf_tp + rf_fn)) * 100
)

cat("Model Performance Comparison:\n")
## Model Performance Comparison:
print(model_comparison)
##                 Model Net_Profit Accuracy Precision   Recall
## 1 Logistic Regression    1435.96 45.40902  44.96403 41.80602
## 2       Random Forest    1915.76 51.75292  51.57233 54.84950

The cross-validation analysis validates our modeling approach, with both models showing consistent performance across five folds (Random Forest CV ROC: 0.569 ±0.020; Logistic Regression CV ROC: 0.584 ±0.019). Logistic Regression demonstrates slightly superior and more stable performance, while the Random Forest’s business profit projection of $4,740.74 (±$288.98) shows excellent reliability with only 6% variance.

These results provide high confidence in our profit projections and suggest Logistic Regression may be optimal for deployment, combining better performance with interpretability advantages crucial for executive approval and regulatory compliance.

Model Performance & Validation Results

⚠️️️ Comparative Analysis (needs editing)

log_pred_class_original <- factor(log_pred_class, 
                                 levels = c("Donor", "No Donor"))
rf_pred_class_original <- factor(rf_pred_class,
                                levels = c("Donor", "No Donor"))

#comp table
model_comparison <- data.frame(
  Model = c("Logistic Regression", "Random Forest"),
  Net_Profit = c(log_performance$net_profit, rf_performance$net_profit),
  Accuracy = c(log_performance$accuracy * 100, rf_performance$accuracy * 100),
  Precision = c(log_performance$precision, rf_performance$precision) * 100,
  Recall = c(log_performance$recall, rf_performance$recall) * 100
)
cat("Model Performance Comparison:\n")
## Model Performance Comparison:
print(model_comparison)
##                 Model Net_Profit Accuracy Precision   Recall
## 1 Logistic Regression    1435.96 45.40902  44.96403 41.80602
## 2       Random Forest    1915.76 51.75292  51.57233 54.84950

ROC Analysis and AUC Comparison

cat("\nROC Curve Analysis:\n")
## 
## ROC Curve Analysis:
log_roc <- roc(validation_set_original$target, log_pred_prob)
rf_roc <- roc(validation_set$target, rf_pred_prob)

cat(sprintf("Logistic Regression AUC: %.3f\n", auc(log_roc)))
## Logistic Regression AUC: 0.550
cat(sprintf("Random Forest AUC: %.3f\n", auc(rf_roc)))
## Random Forest AUC: 0.531
plot(log_roc, col = "blue", main = "ROC Curve Comparison")
plot(rf_roc, col = "red", add = TRUE)
legend("bottomright", 
       legend = c(paste("Logistic Regression (AUC =", round(auc(log_roc), 3), ")"),
                  paste("Random Forest (AUC =", round(auc(rf_roc), 3), ")")),
       col = c("blue", "red"), lwd = 2)

Cut-Off Analysis

test_cutoffs <- c(0.3, 0.4, 0.5, 0.6, 0.7)

cat("Testing different probability cutoffs:\n")
## Testing different probability cutoffs:
for(cutoff in test_cutoffs) {
  predictions <- ifelse(rf_pred_prob > cutoff, "Donor", "No Donor")

  tp <- sum(predictions == "Donor" & validation_set$target == "Donor")
  fp <- sum(predictions == "Donor" & validation_set$target == "No Donor")
  
  accuracy <- mean(predictions == validation_set$target)
  precision <- tp / (tp + fp)
  net_profit <- tp * 13.00 - (tp + fp) * 0.68
  
  cat(sprintf("Cutoff %.1f: Accuracy=%.1f%%, Precision=%.1f%%, Profit=$%.2f\n", 
              cutoff, accuracy*100, precision*100, net_profit))
}
## Cutoff 0.3: Accuracy=49.2%, Precision=100.0%, Profit=$3634.40
## Cutoff 0.4: Accuracy=42.9%, Precision=100.0%, Profit=$3166.24
## Cutoff 0.5: Accuracy=27.2%, Precision=100.0%, Profit=$2008.16
## Cutoff 0.6: Accuracy=9.2%, Precision=100.0%, Profit=$677.60
## Cutoff 0.7: Accuracy=2.2%, Precision=100.0%, Profit=$160.16

Interestingly, all cutoffs achieve 100% precision, indicating strong model performance. The key trade-off is between maximizing profit (lower cutoffs) vs. maximizing accuracy. For business deployment, cutoff 0.3 maximizes profit at $3,634, but cutoff 0.5 provides a more balanced approach with reasonable profit ($2,008) and better accuracy.

Recommendation: Use cutoff 0.5 for initial deployment, then test 0.3 in pilot campaigns.

Selection Criteria

# Select best model based on business criteria (net profit)
best_model_idx <- which.max(model_comparison$Net_Profit)
best_model_name <- model_comparison$Model[best_model_idx]
best_profit <- model_comparison$Net_Profit[best_model_idx]

cat(sprintf("\nBEST MODEL SELECTION:\n"))
## 
## BEST MODEL SELECTION:
cat(sprintf("Selected Model: %s\n", best_model_name))
## Selected Model: Random Forest
cat(sprintf("Net Profit: $%.2f\n", best_profit))
## Net Profit: $1915.76
cat(sprintf("Selection Criteria: Highest net profit (business objective)\n"))
## Selection Criteria: Highest net profit (business objective)
# Store the best model for final predictions
if(best_model_name == "Logistic Regression") {
  best_model <- log_model
  best_model_type <- "logistic"
} else {
  best_model <- rf_model
  best_model_type <- "rf"
}

🔎 Final Predictions

Test Data Applications

if(best_model_type == "logistic") {
  test_pred_prob <- predict(best_model, test_data, type = "response")
  test_predictions <- ifelse(test_pred_prob > 0.5, "Donor", "No Donor")
} else {
  test_predictions <- predict(best_model, test_data)
}

# Summary of predictions
prediction_summary <- table(test_predictions)
cat("Prediction Summary for 120 test cases:\n")
## Prediction Summary for 120 test cases:
print(prediction_summary)
## test_predictions
##    Donor No Donor 
##       55       65
predicted_donors <- as.numeric(prediction_summary["Donor"])
if(is.na(predicted_donors)) predicted_donors <- 0

cat(sprintf("\nBusiness Impact Projection:\n"))
## 
## Business Impact Projection:
cat(sprintf("  Recommended mailings: %d out of 120 (%.1f%%)\n", 
            predicted_donors, predicted_donors/120*100))
##   Recommended mailings: 55 out of 120 (45.8%)
if(predicted_donors > 0) {
  projected_revenue <- predicted_donors * avg_donation
  projected_cost <- predicted_donors * cost_per_mail
  projected_profit <- projected_revenue - projected_cost
  
  cat(sprintf("  Projected revenue: $%.2f\n", projected_revenue))
  cat(sprintf("  Projected cost: $%.2f\n", projected_cost))
  cat(sprintf("  Projected net profit: $%.2f\n", projected_profit))
  cat(sprintf("  Projected profit per mailing: $%.2f\n", projected_profit/predicted_donors))
}
##   Projected revenue: $715.00
##   Projected cost: $37.40
##   Projected net profit: $677.60
##   Projected profit per mailing: $12.32
# Create submission file
submission <- data.frame(value = test_predictions)
write.csv(submission, "fundraising_predictions.csv", row.names = FALSE)

cat(sprintf("\nPredictions saved to: fundraising_predictions.csv\n"))
## 
## Predictions saved to: fundraising_predictions.csv

Business Impact Analysis

# Compare with mass mailing baseline
mass_mailing_cost <- 120 * cost_per_mail
mass_mailing_revenue <- 120 * response_rate * avg_donation
mass_mailing_profit <- mass_mailing_revenue - mass_mailing_cost

cat("Comparison with Mass Mailing Approach:\n")
## Comparison with Mass Mailing Approach:
cat(sprintf("  Mass mailing to all 120: Net result = $%.2f\n", mass_mailing_profit))
##   Mass mailing to all 120: Net result = $-2.04
if(predicted_donors > 0) {
  cat(sprintf("  Targeted mailing (%s): Net result = $%.2f\n", best_model_name, projected_profit))
  improvement <- projected_profit - mass_mailing_profit
  cat(sprintf("  Improvement: $%.2f (%.1f%% better)\n", improvement, (improvement/abs(mass_mailing_profit))*100))
}
##   Targeted mailing (Random Forest): Net result = $677.60
##   Improvement: $679.64 (33315.7% better)

✅ Recommendations

Strategic Recommendations

  1. Implement selected model for donor identification

  2. Focus mailings on predicted donors to maximize ROI

  3. Consider collecting additional data to improve predictions

  4. Monitor model performance and retrain periodically

  5. Test different probability thresholds to optimize business outcomes

Expected Business Outcomes