Adding SCAM with confidence interval as #5 asks

adamesalles · Nov 10, 2024 · fe66e9d · fe66e9d
1 parent f1ec39a
commit fe66e9d
Showing 96 changed files with 460,058 additions and 250,033 deletions.
diff --git a/experiments/adaptive_new_HC/GP_20_nobs.png b/experiments/adaptive_new_HC/GP_20_nobs.png
diff --git a/experiments/adaptive_new_HC/GP_20_nobs_.csv b/experiments/adaptive_new_HC/GP_20_nobs_.csv
diff --git a/experiments/adaptive_new_HC/SCAM_10_nobs.png b/experiments/adaptive_new_HC/SCAM_10_nobs.png
diff --git a/experiments/adaptive_new_HC/SCAM_10_nobs_.csv b/experiments/adaptive_new_HC/SCAM_10_nobs_.csv
diff --git a/experiments/adaptive_new_HC/SCAM_12_nobs.png b/experiments/adaptive_new_HC/SCAM_12_nobs.png
diff --git a/experiments/adaptive_new_HC/SCAM_12_nobs_.csv b/experiments/adaptive_new_HC/SCAM_12_nobs_.csv
diff --git a/experiments/adaptive_new_HC/SCAM_14_nobs.png b/experiments/adaptive_new_HC/SCAM_14_nobs.png
diff --git a/experiments/adaptive_new_HC/SCAM_14_nobs_.csv b/experiments/adaptive_new_HC/SCAM_14_nobs_.csv
diff --git a/experiments/adaptive_new_HC/SCAM_16_nobs.png b/experiments/adaptive_new_HC/SCAM_16_nobs.png
diff --git a/experiments/adaptive_new_HC/SCAM_16_nobs_.csv b/experiments/adaptive_new_HC/SCAM_16_nobs_.csv
diff --git a/experiments/adaptive_new_HC/SCAM_18_nobs.png b/experiments/adaptive_new_HC/SCAM_18_nobs.png
diff --git a/experiments/adaptive_new_HC/SCAM_18_nobs_.csv b/experiments/adaptive_new_HC/SCAM_18_nobs_.csv
diff --git a/experiments/adaptive_new_HC/SCAM_20_nobs.png b/experiments/adaptive_new_HC/SCAM_20_nobs.png
diff --git a/experiments/adaptive_new_HC/SCAM_20_nobs_.csv b/experiments/adaptive_new_HC/SCAM_20_nobs_.csv
diff --git a/experiments/adaptive_new_HC/SCGP_20_nobs.png b/experiments/adaptive_new_HC/SCGP_20_nobs.png
diff --git a/experiments/adaptive_new_HC/SCGP_20_nobs_.csv b/experiments/adaptive_new_HC/SCGP_20_nobs_.csv
diff --git a/experiments/adaptive_new_MC/GP_20_nobs.png b/experiments/adaptive_new_MC/GP_20_nobs.png
diff --git a/experiments/adaptive_new_MC/GP_20_nobs_.csv b/experiments/adaptive_new_MC/GP_20_nobs_.csv
diff --git a/experiments/adaptive_new_MC/SCAM_10_nobs.png b/experiments/adaptive_new_MC/SCAM_10_nobs.png
diff --git a/experiments/adaptive_new_MC/SCAM_10_nobs_.csv b/experiments/adaptive_new_MC/SCAM_10_nobs_.csv
diff --git a/experiments/adaptive_new_MC/SCAM_12_nobs.png b/experiments/adaptive_new_MC/SCAM_12_nobs.png
diff --git a/experiments/adaptive_new_MC/SCAM_12_nobs_.csv b/experiments/adaptive_new_MC/SCAM_12_nobs_.csv
diff --git a/experiments/adaptive_new_MC/SCAM_14_nobs.png b/experiments/adaptive_new_MC/SCAM_14_nobs.png
diff --git a/experiments/adaptive_new_MC/SCAM_14_nobs_.csv b/experiments/adaptive_new_MC/SCAM_14_nobs_.csv
diff --git a/experiments/adaptive_new_MC/SCAM_16_nobs.png b/experiments/adaptive_new_MC/SCAM_16_nobs.png
diff --git a/experiments/adaptive_new_MC/SCAM_16_nobs_.csv b/experiments/adaptive_new_MC/SCAM_16_nobs_.csv
diff --git a/experiments/adaptive_new_MC/SCAM_18_nobs.png b/experiments/adaptive_new_MC/SCAM_18_nobs.png
diff --git a/experiments/adaptive_new_MC/SCAM_18_nobs_.csv b/experiments/adaptive_new_MC/SCAM_18_nobs_.csv
diff --git a/experiments/adaptive_new_MC/SCAM_20_nobs.png b/experiments/adaptive_new_MC/SCAM_20_nobs.png
diff --git a/experiments/adaptive_new_MC/SCAM_20_nobs_.csv b/experiments/adaptive_new_MC/SCAM_20_nobs_.csv
diff --git a/experiments/adaptive_new_MC/SCGP_20_nobs.png b/experiments/adaptive_new_MC/SCGP_20_nobs.png
diff --git a/experiments/adaptive_new_MC/SCGP_20_nobs_.csv b/experiments/adaptive_new_MC/SCGP_20_nobs_.csv
diff --git a/experiments/uniform_new_HC/GP_20_nobs.png b/experiments/uniform_new_HC/GP_20_nobs.png
diff --git a/experiments/uniform_new_HC/GP_20_nobs_.csv b/experiments/uniform_new_HC/GP_20_nobs_.csv
diff --git a/experiments/uniform_new_HC/SCAM_10_nobs.png b/experiments/uniform_new_HC/SCAM_10_nobs.png
diff --git a/experiments/uniform_new_HC/SCAM_10_nobs_.csv b/experiments/uniform_new_HC/SCAM_10_nobs_.csv
diff --git a/experiments/uniform_new_HC/SCAM_12_nobs.png b/experiments/uniform_new_HC/SCAM_12_nobs.png
diff --git a/experiments/uniform_new_HC/SCAM_12_nobs_.csv b/experiments/uniform_new_HC/SCAM_12_nobs_.csv
diff --git a/experiments/uniform_new_HC/SCAM_14_nobs.png b/experiments/uniform_new_HC/SCAM_14_nobs.png
diff --git a/experiments/uniform_new_HC/SCAM_14_nobs_.csv b/experiments/uniform_new_HC/SCAM_14_nobs_.csv
diff --git a/experiments/uniform_new_HC/SCAM_16_nobs.png b/experiments/uniform_new_HC/SCAM_16_nobs.png
diff --git a/experiments/uniform_new_HC/SCAM_16_nobs_.csv b/experiments/uniform_new_HC/SCAM_16_nobs_.csv
diff --git a/experiments/uniform_new_HC/SCAM_18_nobs.png b/experiments/uniform_new_HC/SCAM_18_nobs.png
diff --git a/experiments/uniform_new_HC/SCAM_18_nobs_.csv b/experiments/uniform_new_HC/SCAM_18_nobs_.csv
diff --git a/experiments/uniform_new_HC/SCAM_20_nobs.png b/experiments/uniform_new_HC/SCAM_20_nobs.png
diff --git a/experiments/uniform_new_HC/SCAM_20_nobs_.csv b/experiments/uniform_new_HC/SCAM_20_nobs_.csv
diff --git a/experiments/uniform_new_HC/SCGP_20_nobs.png b/experiments/uniform_new_HC/SCGP_20_nobs.png
diff --git a/experiments/uniform_new_HC/SCGP_20_nobs_.csv b/experiments/uniform_new_HC/SCGP_20_nobs_.csv
diff --git a/experiments/uniform_new_MC/GP_10_nobs.png b/experiments/uniform_new_MC/GP_10_nobs.png
diff --git a/experiments/uniform_new_MC/GP_10_nobs_.csv b/experiments/uniform_new_MC/GP_10_nobs_.csv
diff --git a/experiments/uniform_new_MC/GP_12_nobs.png b/experiments/uniform_new_MC/GP_12_nobs.png
diff --git a/experiments/uniform_new_MC/GP_12_nobs_.csv b/experiments/uniform_new_MC/GP_12_nobs_.csv
diff --git a/experiments/uniform_new_MC/GP_14_nobs.png b/experiments/uniform_new_MC/GP_14_nobs.png
diff --git a/experiments/uniform_new_MC/GP_14_nobs_.csv b/experiments/uniform_new_MC/GP_14_nobs_.csv
diff --git a/experiments/uniform_new_MC/GP_16_nobs.png b/experiments/uniform_new_MC/GP_16_nobs.png
diff --git a/experiments/uniform_new_MC/GP_16_nobs_.csv b/experiments/uniform_new_MC/GP_16_nobs_.csv
diff --git a/experiments/uniform_new_MC/GP_18_nobs.png b/experiments/uniform_new_MC/GP_18_nobs.png
diff --git a/experiments/uniform_new_MC/GP_18_nobs_.csv b/experiments/uniform_new_MC/GP_18_nobs_.csv
diff --git a/experiments/uniform_new_MC/GP_20_nobs.png b/experiments/uniform_new_MC/GP_20_nobs.png
diff --git a/experiments/uniform_new_MC/GP_20_nobs_.csv b/experiments/uniform_new_MC/GP_20_nobs_.csv
diff --git a/experiments/uniform_new_MC/GP_6_nobs.png b/experiments/uniform_new_MC/GP_6_nobs.png
diff --git a/experiments/uniform_new_MC/GP_6_nobs_.csv b/experiments/uniform_new_MC/GP_6_nobs_.csv
diff --git a/experiments/uniform_new_MC/GP_8_nobs.png b/experiments/uniform_new_MC/GP_8_nobs.png
diff --git a/experiments/uniform_new_MC/GP_8_nobs_.csv b/experiments/uniform_new_MC/GP_8_nobs_.csv
diff --git a/experiments/uniform_new_MC/SCAM_10_nobs.png b/experiments/uniform_new_MC/SCAM_10_nobs.png
diff --git a/experiments/uniform_new_MC/SCAM_10_nobs_.csv b/experiments/uniform_new_MC/SCAM_10_nobs_.csv
diff --git a/experiments/uniform_new_MC/SCAM_12_nobs.png b/experiments/uniform_new_MC/SCAM_12_nobs.png
diff --git a/experiments/uniform_new_MC/SCAM_12_nobs_.csv b/experiments/uniform_new_MC/SCAM_12_nobs_.csv
diff --git a/experiments/uniform_new_MC/SCAM_14_nobs.png b/experiments/uniform_new_MC/SCAM_14_nobs.png
diff --git a/experiments/uniform_new_MC/SCAM_14_nobs_.csv b/experiments/uniform_new_MC/SCAM_14_nobs_.csv
diff --git a/experiments/uniform_new_MC/SCAM_16_nobs.png b/experiments/uniform_new_MC/SCAM_16_nobs.png
diff --git a/experiments/uniform_new_MC/SCAM_16_nobs_.csv b/experiments/uniform_new_MC/SCAM_16_nobs_.csv
diff --git a/experiments/uniform_new_MC/SCAM_18_nobs.png b/experiments/uniform_new_MC/SCAM_18_nobs.png
diff --git a/experiments/uniform_new_MC/SCAM_18_nobs_.csv b/experiments/uniform_new_MC/SCAM_18_nobs_.csv
diff --git a/experiments/uniform_new_MC/SCAM_20_nobs.png b/experiments/uniform_new_MC/SCAM_20_nobs.png
diff --git a/experiments/uniform_new_MC/SCAM_20_nobs_.csv b/experiments/uniform_new_MC/SCAM_20_nobs_.csv
diff --git a/experiments/uniform_new_MC/SCGP_10_nobs.png b/experiments/uniform_new_MC/SCGP_10_nobs.png
diff --git a/experiments/uniform_new_MC/SCGP_10_nobs_.csv b/experiments/uniform_new_MC/SCGP_10_nobs_.csv
diff --git a/experiments/uniform_new_MC/SCGP_12_nobs.png b/experiments/uniform_new_MC/SCGP_12_nobs.png
diff --git a/experiments/uniform_new_MC/SCGP_12_nobs_.csv b/experiments/uniform_new_MC/SCGP_12_nobs_.csv
diff --git a/experiments/uniform_new_MC/SCGP_14_nobs.png b/experiments/uniform_new_MC/SCGP_14_nobs.png
diff --git a/experiments/uniform_new_MC/SCGP_14_nobs_.csv b/experiments/uniform_new_MC/SCGP_14_nobs_.csv
diff --git a/experiments/uniform_new_MC/SCGP_16_nobs.png b/experiments/uniform_new_MC/SCGP_16_nobs.png
diff --git a/experiments/uniform_new_MC/SCGP_16_nobs_.csv b/experiments/uniform_new_MC/SCGP_16_nobs_.csv
diff --git a/experiments/uniform_new_MC/SCGP_18_nobs.png b/experiments/uniform_new_MC/SCGP_18_nobs.png
diff --git a/experiments/uniform_new_MC/SCGP_18_nobs_.csv b/experiments/uniform_new_MC/SCGP_18_nobs_.csv
diff --git a/experiments/uniform_new_MC/SCGP_20_nobs.png b/experiments/uniform_new_MC/SCGP_20_nobs.png
diff --git a/experiments/uniform_new_MC/SCGP_20_nobs_.csv b/experiments/uniform_new_MC/SCGP_20_nobs_.csv
diff --git a/experiments/uniform_new_MC/SCGP_6_nobs.png b/experiments/uniform_new_MC/SCGP_6_nobs.png
diff --git a/experiments/uniform_new_MC/SCGP_6_nobs_.csv b/experiments/uniform_new_MC/SCGP_6_nobs_.csv
diff --git a/experiments/uniform_new_MC/SCGP_8_nobs.png b/experiments/uniform_new_MC/SCGP_8_nobs.png
diff --git a/experiments/uniform_new_MC/SCGP_8_nobs_.csv b/experiments/uniform_new_MC/SCGP_8_nobs_.csv
diff --git a/src/experiments/reducing-points-experiment.py b/src/experiments/reducing-points-experiment.py
@@ -34,7 +34,7 @@
 
     scenarios = {
         "nobs": list(range(6, 21, 2)),
-        # "nobs": [10, 20],
+        # "nobs": [20],
     }
 
     mse_df = pd.DataFrame(columns=["model", "scenario", "dataset", "mse"])

diff --git a/src/scam/plot.py b/src/scam/plot.py
@@ -1,8 +1,6 @@
-import pandas as pd
 import numpy as np
 import matplotlib
 import matplotlib.pyplot as plt
-import torch
 import pathlib
 
 # Define colors
@@ -24,14 +22,18 @@
         "adaptive_new_MC": DATA_PATH / "Gaussian_logCA0_adaptive_J=20_MC.csv",
     }
 
+
 # Loading data
-def load_data(path: str) -> tuple:
-    data = np.loadtxt(path, delimiter=",", skiprows=1, dtype=np.float32)
+def load_data(path: str, se: bool = False) -> tuple:
+    data = np.loadtxt(path, delimiter=",", skiprows=2, dtype=np.float32)
     x_train = data[:, 0]
     y_train = data[:, 1]
-
+    if se:
+        se_train = data[:, 2]
+        return x_train, y_train, se_train
     return x_train, y_train
 
+
 def remove_points(train_x, train_y, nobs):
     """
     Remove points from the dataset equally spaced
@@ -51,7 +53,7 @@ def remove_points(train_x, train_y, nobs):
         else:
             x_true, y_true = true_MC_x[1:], true_MC_y[1:]
 
-        for nobs in range(10, 21, 2) + [21]:
+        for nobs in range(10, 21, 2):
             fig, ax = plt.subplots(1, 1, figsize=(5, 4), tight_layout=True)
             name = f"{model}_{nobs}_nobs"
             io_path = PATH / "experiments" / str(dataset_name)
@@ -60,10 +62,15 @@ def remove_points(train_x, train_y, nobs):
             train_x, train_y = load_data(dataset_path)
             train_x, train_y = remove_points(train_x, train_y, nobs)
 
-            test_x, test_y = load_data(io_path / f'{name}_.csv')
+            test_x, test_y, test_se = load_data(io_path / f'{name}_.csv',
+                                                se=True)
             mse_loss = np.mean((test_y - y_true)**2)
 
             ax.plot(x_true, y_true, color=COLORS[0])
+            ax.fill_between(
+                test_x, test_y + (2*test_se), test_y - (2*test_se),
+                color=COLORS[0], alpha=0.5
+            )
             ax.plot(train_x, train_y, "k*")
             ax.legend(["Observed Values", "Mean", "Confidence"])
             ax.set_title("Function values")
@@ -77,10 +84,11 @@ def remove_points(train_x, train_y, nobs):
                     bbox=dict(facecolor=COLORS[0], alpha=0.5))
 
             fig.suptitle("Shape Constrained Additive Model (SCAM)",
-               fontsize=18,
-               color=COLORS[0],
-               fontweight='bold')
+                         fontsize=18,
+                         color=COLORS[0],
+                         fontweight='bold')
 
             # plt.show()
-            fig.savefig(io_path / (name + ".png"), dpi=600, bbox_inches="tight")
+            fig.savefig(io_path / (name + ".png"), dpi=600,
+                        bbox_inches="tight")
             plt.close(fig)
diff --git a/src/scam/scam.R b/src/scam/scam.R
@@ -17,7 +17,7 @@ data_paths <- list(
 )
 
 # Define nobs scenarios
-nobs_scenarios <- append(seq(6, 20, by = 2), 21)
+nobs_scenarios <- seq(6, 20, by = 2)
 
 # Data frame to store results
 results <- data.frame(model = character(), scenario = integer(), dataset = character(), mse = numeric())
@@ -28,7 +28,7 @@ for (dataset_name in names(data_paths)) {
   dataset_path <- data_paths[[dataset_name]]
 
   # Read data
-  data <- read_csv(dataset_path)
+  data <- read_csv(dataset_path) %>% filter(!row_number() == 1)
   x <- data$a0
   y <- data$lc_a0
 
@@ -54,23 +54,17 @@ for (dataset_name in names(data_paths)) {
     if (!is.null(model)) {
       # Generate predictions on a regular grid from 0 to 1
       grid <- data.frame(x_sampled = x_true)
-      grid$y_pred <- predict(model, newdata = grid)
+      p <- predict(model, newdata = grid, se.fit = TRUE)
+      grid$y_pred <- p$fit
+      grid$se <- p$se.fit
+
+      plot(model)
       # print("Got here")
 
       # Save predictions
       output_file <- paste0("../../experiments/", dataset_name, "/SCAM_", nobs , "_nobs_",  ".csv")
       write.csv(grid, output_file, row.names = FALSE)
 
-      # Calculate Mean Squared Error (MSE) on the sampled data
-      #y_pred_sampled <- predict(model, newdata = data.frame(x = x_sampled))
-      #mse <- mean((y_sampled - y_pred_sampled)^2)
-
-      # Store results
-      #results <- rbind(
-      #  results,
-      #  data.frame(model = "SCAM", scenario = nobs, dataset = dataset_name, mse = mse)
-      #)
-
       cat("Finished SCAM model on", dataset_name, "with nobs =", nobs, "\n")
       rm(grid)
     }

diff --git a/src/scgp/utils.py b/src/scgp/utils.py
@@ -3,6 +3,7 @@
 import matplotlib
 import pathlib
 import numpy as np
+import pandas as pd
 import matplotlib.pyplot as plt
 from tqdm import tqdm, trange
 from scgp.models import GP, SCGP
@@ -25,7 +26,7 @@
 
 # Loading data
 def load_data(path: str) -> tuple:
-    data = np.loadtxt(path, delimiter=",", skiprows=1, dtype=np.float32)
+    data = np.loadtxt(path, delimiter=",", skiprows=2, dtype=np.float32)
     x_train = torch.from_numpy(data[:, 0])
     y_train = torch.stack(
         [torch.from_numpy(data[:, 1]), torch.from_numpy(data[:, 2])], -1
@@ -154,15 +155,17 @@ def save_plot_scgp(
         color=COLORS[1]
     )
     y_ax.plot(train_x.numpy(), train_y[:, 0].numpy(), "k*")
-    y_ax.legend(["Observed Values", "Mean", "Confidence"])
+
+    y_ax.legend(["Observed Values", "Mean", "Confidence"],
+                loc="upper left")
     y_ax.set_title("Function values")
     y_ax.set_xlim([0, 1])
     # y_ax.set_ylim([-7.5, 12.5])
     y_ax.set_xlabel(r"$\alpha$")
     y_ax.set_ylabel(r"$f_{\boldsymbol{z}}(\alpha)$")
     if plot_mse:
         # show mse in bottom right (larger font)
-        y_ax.text(0.7, 0.05, f"MSE: {mse_loss:.2f}",
+        y_ax.text(0.5, 0.05, f"MSE: {mse_loss:.2f}",
                   transform=y_ax.transAxes,
                     fontsize=16,
                     bbox=dict(facecolor=COLORS[1], alpha=0.5))
@@ -175,15 +178,17 @@ def save_plot_scgp(
     )
     y_prime_ax.plot(train_x.numpy(), train_y[:, 1].numpy(),
                     "k*")
-    y_prime_ax.legend(["Observed Derivatives", "Mean", "Confidence"])
+
+    y_prime_ax.legend(["Observed Derivatives", "Mean", "Confidence"],
+                      loc="upper left")
     y_prime_ax.set_title(r"Derivatives with respect to $\alpha$")
     y_prime_ax.set_xlim([0, 1])
 
     y_prime_ax.set_xlabel(r"$\alpha$")
     y_prime_ax.set_ylabel(r"$\frac{\mathrm{d}}{\mathrm{d}\alpha}f_{\boldsymbol{z}}(\alpha)$")
     if plot_mse:
         # show mse
-        y_prime_ax.text(0.7, 0.05, f"MSE: {mse_loss_prime:.2f}",
+        y_prime_ax.text(0.5, 0.05, f"MSE: {mse_loss_prime:.2f}",
                   transform=y_prime_ax.transAxes,
                     fontsize=16,
                     bbox=dict(facecolor=COLORS[1], alpha=0.5))
@@ -201,7 +206,21 @@ def save_plot_scgp(
     # f.patch.set_linewidth(6)
     # f.patch.set_edgecolor(COLORS[1])
 
-
+    # save csv with test data
+    if save_path:
+        pd.DataFrame(
+            {
+                "x": test_x.numpy(),
+                "f": test_y[:, 0].numpy(),
+                "f_prime": test_y[:, 1].numpy(),
+                "mean": mean[:, 0].numpy(),
+                "mean_prime": mean[:, 1].numpy(),
+                "lower": lower[:, 0].numpy(),
+                "lower_prime": lower[:, 1].numpy(),
+                "upper": upper[:, 0].numpy(),
+                "upper_prime": upper[:, 1].numpy(),
+            }
+        ).to_csv(save_path / (name + "_.csv"), index=False)
     f.savefig(save_path / (name + ".png"), dpi=600, bbox_inches="tight")
     plt.close(f)
 
@@ -305,7 +324,7 @@ def save_plot_gp(
     y_ax.set_ylabel(r"$f_{\boldsymbol{z}}(\alpha)$")
     if plot_mse:
         # show mse
-        y_ax.text(0.7, 0.05, f"MSE: {mse_loss:.2f}",
+        y_ax.text(0.5, 0.05, f"MSE: {mse_loss:.2f}",
                   transform=y_ax.transAxes,
                   fontsize=16,
                   bbox=dict(facecolor=COLORS[0], alpha=0.5))
@@ -321,7 +340,16 @@ def save_plot_gp(
 
     # f.patch.set_linewidth(6)
     # f.patch.set_edgecolor(COLORS[0])
-
+    if save_path:
+        pd.DataFrame(
+            {
+                "x": test_x.numpy(),
+                "f": test_y.numpy(),
+                "mean": mean.numpy(),
+                "lower": lower.numpy(),
+                "upper": upper.numpy(),
+            }
+        ).to_csv(save_path / (name + "_.csv"), index=False)
     f.savefig(save_path / (name + ".png"), dpi=600, bbox_inches="tight")
     plt.close(f)