rm(list = ls())
library(MASS)
#[Link]("robustbase")
library(robustbase)
library(ggplot2)
# Set seed for reproducibility
[Link](811)
# Generate independent variable X (100 observations)
X <- runif(1000, min = 30, max = 45)
X1 <- runif(1000, min = 30, max = 45)
# Generate dependent variable Y with noise (related to X)
Y <- 3*X + 2 + rnorm(1000, mean = 0, sd = 2)
Y1 <- 3*X + 2 + rnorm(1000, mean = 0, sd = 2)
# Introduce outliers in X and Y
X[615] <- 400
X[510] <- 90
X[230] <- 450
X[110] <- 105
X[21] <- 355
# Introduce outliers in X
data1 = [Link](X=X, Y=Y)
data1
Y1[79] <- 366
Y1[200] <- 375
Y1[251] <- 450
Y1[399] <- 390
Y1[402] <- 460
# Corresponding outliers in Y
# Create a dataframe named 'data' with X and Y
data2 <- [Link](X = X1, Y = Y1)
# Display the dataframe
print(data2)
# Get the directory path of the R script
# Using psych package for descriptive statistics
#[Link]("psych") # Install if not already installed
library(psych)
# Descriptive statistics with psych package for data1
des1 = describe(data1)
des1
# Descriptive statistics with psych package for data2
des2 = describe(data2)
des2
# Load ggplot2 library
library(ggplot2)
# Create a box plot for data1
ggplot(data1, aes(x = X, y = "")) +
geom_boxplot(fill = "8") + # Customize box plot appearance
theme_minimal() + # Customize plot theme (optional)
labs(y = "Y", x = "X (Data 1)") # Customize axis labels (optional)
# Create a box plot for data1
ggplot(data2, aes(x = "", y = Y)) +
geom_boxplot(fill = "8") + # Customize box plot appearance
theme_minimal() + # Customize plot theme (optional)
labs(x = "X", y = "Y (Data 2)") # Customize axis labels (optional)
### Linear
ols_model = lm(data1$Y~data1$X, data = data1)
summary(ols_model)
m_est_model = rlm(data1$Y~data1$X, data = data1)
summary(m_est_model)
# Fit Least Median of Squares (LMS) regression
lms_model <- lqs(Y ~ X, data = data1)
summary(lms_model)
# Fit Least Absolute Deviation (LAD) regression (robust to outliers)
lad_model <- lmrob(Y ~ X, data = data1)
summary(lad_model)
## LTS
lts_model = robustbase::ltsReg(data1$Y~data1$X, data = data1)
summary(lts_model)
# MM
mm_model = robustbase::lmrob(data1$Y~data1$X, data = data1)
summary(mm_model)
# Fit Huber Regression
huber_model <- lmrob(data1$Y~data1$X, data = data1, method = "MM") # Use method =
"MM" for Huber Regression
summary(huber_model)
# Fit RANSAC Regression
ransac_model <- rlm(data1$Y~data1$X, data = data1, method = "MM") # Use method =
"MM" for RANSAC Regression
summary(ransac_model)
########
## Linear
ols_model1 = lm(data2$Y~data2$X, data = data2)
summary(ols_model1)
## M estimation
m_est_model1 = rlm(data2$Y~data2$X, data = data2)
summary(m_est_model1)
# Fit Least Median of Squares (LMS) regression
lms_model1 <- lqs(data2$Y~data2$X, data = data2)
# Fit Least Absolute Deviation (LAD) regression (robust to outliers)
lad_model1 <- lmrob(data2$Y~data2$X, data = data2)
## LTS
lts_model1 = robustbase::ltsReg(data2$Y~data2$X, data = data2)
summary(lts_model1)
# MM
mm_model1 = robustbase::lmrob(data2$Y~data2$X, data = data2)
summary(mm_model1)
# Fit Huber Regression
huber_model1 <- lmrob(data2$Y~data2$X, data = data2, method = "MM") # Use method =
"MM" for Huber Regression
# Fit RANSAC Regression
ransac_model1 <- rlm(data2$Y~data2$X, data = data2, method = "MM") # Use method =
"MM" for RANSAC Regression
###########
# Plot the dataset and regression lines
plot(data1$X, data1$Y, main = "Linear Model (Data 1)", xlab = "X", ylab = "Y",
pch = 1)
# Add linear regression line (OLS)
abline(ols_model, col = "blue", lwd = 2,lty=1)
# Plot the dataset and regression lines
plot(data2$X, data2$Y, main = "Linear Model (Data 2)", xlab = "X", ylab = "Y",
pch = 1)
# Add linear regression line (OLS)
abline(ols_model1, col = "blue", lwd = 2,lty=1)
########
# Plot the dataset and regression lines
plot(data1$X, data1$Y, main = "LMS Model (Data 1)", xlab = "X", ylab = "Y", pch =
1)
abline(lms_model, col = "red", lwd = 2, lty = 1)
# Plot the dataset and regression lines
plot(data2$X, data2$Y, main = "LMS Model (Data 2)", xlab = "X", ylab = "Y", pch =
1)
abline(lms_model1, col = "red", lwd = 2, lty = 1)
###########
# Plot the dataset and regression lines
plot(data1$X, data1$Y, main = "LAD Model (Data 1)", xlab = "X", ylab = "Y", pch =
1)
abline(lad_model, col = "cyan", lwd = 2, lty = 1)
# Plot the dataset and regression lines
plot(data2$X, data2$Y, main = "LAD Model (Data 2)", xlab = "X", ylab = "Y", pch =
1)
abline(lad_model1, col = "cyan", lwd = 2, lty = 1)
###########
# Plot the dataset and regression lines
plot(data1$X, data1$Y, main = "M Estimatoin (Data 1)", xlab = "X", ylab = "Y",
pch = 1)
# Add M-estimation line
abline(m_est_model, col = "green", lwd = 2, lty = 1)
# Plot the dataset and regression lines
plot(data2$X, data2$Y, main = "M Estimation (Data 2)", xlab = "X", ylab = "Y",
pch = 1)
# Add M-estimation line
abline(m_est_model1, col = "green", lwd = 2, lty = 1)
####################
# Plot the dataset and regression lines
plot(data1$X, data1$Y, main = "LTS Model (Data 1)", xlab = "X", ylab = "Y", pch =
1)
# Add LTS line
abline(lts_model, col = "purple", lwd = 2, lty = 1)
# Plot the dataset and regression lines
plot(data2$X, data2$Y, main = "LTS Model (Data 2)", xlab = "X", ylab = "Y", pch =
1)
# Add LTS line
abline(lts_model1, col = "purple", lwd = 2, lty = 1)
###############
# Plot the dataset and regression lines
plot(data1$X, data1$Y, main = "MM Model (Data 1)", xlab = "X", ylab = "Y", pch =
1)
# Add MM line
abline(mm_model, col = "orange", lwd = 2, lty = 1)
# Plot the dataset and regression lines
plot(data2$X, data2$Y, main = "MM Model (Data 2)", xlab = "X", ylab = "Y", pch =
1)
# Add MM line
abline(mm_model1, col = "orange", lwd = 2, lty = 1)
#################
# Plot the dataset and regression lines
plot(data1$X, data1$Y, main = "Huber Model (Data 1)", xlab = "X", ylab = "Y", pch
= 1)
abline(huber_model, col = "purple", lwd = 2, lty = 1)
# Plot the dataset and regression lines
plot(data2$X, data2$Y, main = "Huber Model (Data 2)", xlab = "X", ylab = "Y", pch
= 1)
abline(huber_model1, col = "purple", lwd = 2, lty = 1)
#################
# Plot the dataset and regression lines
plot(data1$X, data1$Y, main = "RANSAC Model (Data 1)", xlab = "X", ylab = "Y",
pch = 1)
abline(ransac_model, col = "pink", lwd = 2, lty = 1)
# Plot the dataset and regression lines
plot(data2$X, data2$Y, main = "RANSAC Model (Data 2)", xlab = "X", ylab = "Y",
pch = 1)
abline(ransac_model1, col = "pink", lwd = 2, lty = 1)
################
# Plot the dataset and regression lines
plot(data1$X, data1$Y, main = "Robust Regression Models (Data 1)", xlab = "X",
ylab = "Y", pch = 1)
# Add linear regression line (OLS)
abline(ols_model, col = "blue", lwd = 2,lty=1)
abline(lms_model, col = "red", lwd = 2, lty = 2)
abline(lad_model, col = "cyan", lwd = 2, lty = 3)
# Add M-estimation line
abline(m_est_model, col = "green", lwd = 2, lty = 4)
# Add LTS line
abline(lts_model, col = "green", lwd = 2, lty = 5)
# Add MM line
abline(mm_model, col = "orange", lwd = 2, lty = 6)
abline(huber_model, col = "purple", lwd = 2, lty = 7)
abline(ransac_model, col = "pink", lwd = 2, lty = 10)
# Add legend for regression lines
legend("topright", legend = c("OLS", "M-Estimation","LMS","LAD", "LTS",
"MM","Huber","Ransac"),
col = c("blue", "green","red","cyan", "green", "orange","purple","pink"),
lty = c(1, 2, 3, 4),
lwd = 2,
title = "Regression Models",cex=.45)
#################### for data2
# Plot the dataset and regression lines
plot(data2$X, data2$Y, main = "Robust Regression Models (Data 2)", xlab = "X",
ylab = "Y", pch = 1)
abline(ols_model1, col = "red", lwd = 2,lty=1)
abline(m_est_model1, col = "blue", lwd = 2, lty = 2)
abline(lms_model1, col = "gold", lwd = 2, lty = 3)
abline(lad_model1, col = "cyan", lwd = 2, lty = 4)
abline(lts_model1, col = "green", lwd = 2, lty = 5)
abline(mm_model1, col = "orange", lwd = 2, lty = 6)
abline(huber_model1, col = "purple", lwd = 2, lty = 7)
abline(ransac_model1, col = "pink", lwd = 2, lty = 10)
# Add legend for regression lines
legend("topright", legend = c("OLS", "M-Estimation","LMS","LAD", "LTS",
"MM","Huber","Ransac"),
col = c("red", "blue","gold","cyan", "green", "orange","purple","pink"),
lty = c(1, 2, 3, 4),
lwd = 2,
title = "Regression Models",cex=.45)
# Load necessary libraries
library(Metrics) # For calculating MAE and RMSE
# Assuming models are already fitted or defined correctly, replace the placeholders
with your model objects
# Define a function to calculate MAE and RMSE
calculate_errors <- function(actual, predicted) {
mae <- mae(actual, predicted) # Calculate Mean Absolute Error (MAE)
rmse <- rmse(actual, predicted) # Calculate Root Mean Squared Error (RMSE)
return(list(MAE = mae, RMSE = rmse))
}
# Calculate MAE and RMSE for each model
mae_rmse_results <- list()
# Assuming ols_model, m_est_model, lms_model, lad_model, lts_model, mm_model are
correctly defined or fitted
# Replace with your actual model objects
# Ordinary Least Squares (OLS) Linear Regression
ols_predicted <- predict(ols_model, newdata = data1)
mae_rmse_results[["OLS"]] <- calculate_errors(data1$Y, ols_predicted)
# Robust Regression - M-Estimation
m_est_predicted <- predict(m_est_model, newdata = data1)
mae_rmse_results[["M Estimation"]] <- calculate_errors(data1$Y, m_est_predicted)
# Robust Regression - Least Median of Squares (LMS)
lms_predicted <- fitted(lms_model)
mae_rmse_results[["LMS"]] <- calculate_errors(data1$Y, lms_predicted)
# Robust Regression - Least Absolute Deviation (LAD)
lad_predicted <- predict(lad_model, newdata = data1)
mae_rmse_results[["LAD"]] <- calculate_errors(data1$Y, lad_predicted)
# Robust Regression - MM (M-Estimator)
mm_predicted <- fitted(mm_model)
mae_rmse_results[["MM"]] <- calculate_errors(data1$Y, mm_predicted)
# Combine results into a data frame
mae_rmse_df <- [Link](rbind, mae_rmse_results)
# Add model names as row names
rownames(mae_rmse_df) <- names(mae_rmse_results)
# Print the MAE and RMSE for each model
print(mae_rmse_df)
# Load necessary libraries
library(Metrics) # For calculating MAE and RMSE
library(robustbase) # For robust regression (lts)
# Define a function to calculate MAE and RMSE
calculate_errors <- function(actual, predicted) {
mae <- mae(actual, predicted) # Calculate Mean Absolute Error (MAE)
rmse <- rmse(actual, predicted) # Calculate Root Mean Squared Error (RMSE)
return(list(MAE = mae, RMSE = rmse))
}
# Example: Fit an LTS (Least Trimmed Squares) model
lts_model <- robustbase::ltsReg(Y ~ X, data = data1)
# Get fitted values (predictions) from the LTS model
lts_predicted <- fitted(lts_model)
# Calculate MAE and RMSE for the LTS model
mae_rmse_results <- calculate_errors(data1$Y, lts_predicted)
# Print the MAE and RMSE for the LTS model
print(mae_rmse_results)
###################
###################
##################
##################
# Load necessary libraries
library(Metrics) # For calculating MAE and RMSE
# Assuming models are already fitted or defined correctly, replace the placeholders
with your model objects
# Define a function to calculate MAE and RMSE
calculate_errors <- function(actual, predicted) {
mae <- mae(actual, predicted) # Calculate Mean Absolute Error (MAE)
rmse <- rmse(actual, predicted) # Calculate Root Mean Squared Error (RMSE)
return(list(MAE = mae, RMSE = rmse))
}
# Calculate MAE and RMSE for each model
mae_rmse_results <- list()
# Assuming ols_model, m_est_model, lms_model, lad_model, lts_model, mm_model are
correctly defined or fitted
# Replace with your actual model objects
# Ordinary Least Squares (OLS) Linear Regression
ols_predicted <- predict(ols_model, newdata = data2)
mae_rmse_results[["OLS"]] <- calculate_errors(data2$Y, ols_predicted)
# Robust Regression - M-Estimation
m_est_predicted <- predict(m_est_model, newdata = data2)
mae_rmse_results[["M Estimation"]] <- calculate_errors(data2$Y, m_est_predicted)
# Robust Regression - Least Median of Squares (LMS)
lms_predicted <- fitted(lms_model)
mae_rmse_results[["LMS"]] <- calculate_errors(data2$Y, lms_predicted)
# Robust Regression - Least Absolute Deviation (LAD)
lad_predicted <- predict(lad_model, newdata = data2)
mae_rmse_results[["LAD"]] <- calculate_errors(data2$Y, lad_predicted)
# Robust Regression - MM (M-Estimator)
mm_predicted <- fitted(mm_model)
mae_rmse_results[["MM"]] <- calculate_errors(data2$Y, mm_predicted)
# Combine results into a data frame
mae_rmse_df <- [Link](rbind, mae_rmse_results)
# Add model names as row names
rownames(mae_rmse_df) <- names(mae_rmse_results)
# Print the MAE and RMSE for each model
print(mae_rmse_df)
# Load necessary libraries
library(Metrics) # For calculating MAE and RMSE
library(robustbase) # For robust regression (lts)
# Define a function to calculate MAE and RMSE
calculate_errors <- function(actual, predicted) {
mae <- mae(actual, predicted) # Calculate Mean Absolute Error (MAE)
rmse <- rmse(actual, predicted) # Calculate Root Mean Squared Error (RMSE)
return(list(MAE = mae, RMSE = rmse))
}
# Example: Fit an LTS (Least Trimmed Squares) model
lts_model <- robustbase::ltsReg(Y ~ X, data = data2)
# Get fitted values (predictions) from the LTS model
lts_predicted <- fitted(lts_model)
# Calculate MAE and RMSE for the LTS model
mae_rmse_results <- calculate_errors(data2$Y, lts_predicted)
# Print the MAE and RMSE for the LTS model
print(mae_rmse_results)
# Fit Ordinary Least Squares (OLS) Linear Regression Model
ols_model_data1 <- lm(Y ~ X, data = data1)
# Get summary of OLS model
summary_data1 <- summary(ols_model_data1)
# Extract R-squared and adjusted R-squared
r_squared_data1 <- summary_data1$[Link]
adj_r_squared_data1 <- summary_data1$[Link]
# Calculate AIC and BIC
aic_data1 <- AIC(ols_model_data1)
bic_data1 <- BIC(ols_model_data1)
# Print results for data1
cat("R-squared (data1):", r_squared_data1, "\n")
cat("Adjusted R-squared (data1):", adj_r_squared_data1, "\n")
cat("AIC (data1):", aic_data1, "\n")
cat("BIC (data1):", bic_data1, "\n")
# Fit Ordinary Least Squares (OLS) Linear Regression Model
ols_model_data2 <- lm(Y ~ X, data = data2)
# Get summary of OLS model
summary_data2 <- summary(ols_model_data2)
# Extract R-squared and adjusted R-squared
r_squared_data2 <- summary_data2$[Link]
adj_r_squared_data2 <- summary_data2$[Link]
# Calculate AIC and BIC
aic_data2 <- AIC(ols_model_data2)
bic_data2 <- BIC(ols_model_data2)
# Print results for data2
cat("R-squared (data2):", r_squared_data2, "\n")
cat("Adjusted R-squared (data2):", adj_r_squared_data2, "\n")
cat("AIC (data2):", aic_data2, "\n")
cat("BIC (data2):", bic_data2, "\n")
