nlpd score added

examples/example_kernels.py

@@ -183,3 +183,57 @@ def deep_multi_task_kernel(x1, x2, hps):  # pragma: no cover
     d = get_distance_matrix(x1_nn, x2_nn)
     k = signal_var * matern_kernel_diff1(d, length_scale)
     return k
+
+
+
+import torch
+import torch.nn as nn
+import torch.optim as optim
+
+
+# Define a simple neural network to warp a 3D space
+class WarpNet(nn.Module):
+    def __init__(self, input_dim=3, hidden_dim=64, output_dim=3):
+        super(WarpNet, self).__init__()
+
+        # Define the architecture
+        self.fc1 = nn.Linear(input_dim, hidden_dim)
+        self.fc2 = nn.Linear(hidden_dim, hidden_dim)
+        self.fc3 = nn.Linear(hidden_dim, output_dim)
+
+        # Activation functions
+        self.relu = nn.ReLU()
+
+    def forward(self, x):
+        # Pass input through the layers
+        x = self.relu(self.fc1(x))  # Input layer to hidden layer 1
+        x = self.relu(self.fc2(x))  # Hidden layer 1 to hidden layer 2
+        x = self.fc3(x)  # Hidden layer 2 to output layer
+        return x
+
+
+# Initialize the network
+model = WarpNet()
+
+# Example usage: Warping a 3D point (e.g., [x, y, z])
+#points = torch.tensor([[1.0, 2.0, 3.0], [4.0, 5.0, 6.0]])  # Example input
+#warped_points = model(points)  # Get warped points
+#print("Warped Points: ", warped_points)
+
+# Loss function and optimizer for training
+#criterion = nn.MSELoss()
+#optimizer = optim.Adam(model.parameters(), lr=0.001)
+
+# Dummy target points (assuming you know where the warped points should go)
+#target_points = torch.tensor([[0.5, 1.5, 2.5], [3.5, 4.5, 5.5]])
+
+# Training loop (single step for simplicity)
+#optimizer.zero_grad()
+#output = model(points)
+#loss = criterion(output, target_points)  # Calculate loss
+#loss.backward()  # Backpropagation
+#optimizer.step()  # Update the weights
+
+#print("Updated Warped Points: ", model(points))
+
+

fvgp/deep_kernel_network.py

@@ -33,3 +33,54 @@ def get_weights(self):
 
     def get_biases(self):
         return self.layer1.bias, self.layer2.bias, self.layer3.bias
+
+
+import torch
+import torch.nn as nn
+import torch.optim as optim
+
+
+# Define a simple neural network to warp a 3D space
+class WarpNet(nn.Module):
+    def __init__(self, input_dim=3, hidden_dim=64, output_dim=3):
+        super(WarpNet, self).__init__()
+
+        # Define the architecture
+        self.fc1 = nn.Linear(input_dim, hidden_dim)
+        self.fc2 = nn.Linear(hidden_dim, hidden_dim)
+        self.fc3 = nn.Linear(hidden_dim, output_dim)
+
+        # Activation functions
+        self.relu = nn.ReLU()
+
+    def forward(self, x):
+        # Pass input through the layers
+        x = self.relu(self.fc1(x))  # Input layer to hidden layer 1
+        x = self.relu(self.fc2(x))  # Hidden layer 1 to hidden layer 2
+        x = self.fc3(x)  # Hidden layer 2 to output layer
+        return x
+
+
+# Initialize the network
+model = WarpNet()
+
+# Example usage: Warping a 3D point (e.g., [x, y, z])
+#points = torch.tensor([[1.0, 2.0, 3.0], [4.0, 5.0, 6.0]])  # Example input
+#warped_points = model(points)  # Get warped points
+#print("Warped Points: ", warped_points)
+
+# Loss function and optimizer for training
+#criterion = nn.MSELoss()
+#optimizer = optim.Adam(model.parameters(), lr=0.001)
+
+# Dummy target points (assuming you know where the warped points should go)
+#target_points = torch.tensor([[0.5, 1.5, 2.5], [3.5, 4.5, 5.5]])
+
+# Training loop (single step for simplicity)
+#optimizer.zero_grad()
+#output = model(points)
+#loss = criterion(output, target_points)  # Calculate loss
+#loss.backward()  # Backpropagation
+#optimizer.step()  # Update the weights
+
+#print("Updated Warped Points: ", model(points))

fvgp/gp.py

@@ -1301,7 +1301,7 @@ def crps(self, x_test, y_test):
 
         mean = self.posterior_mean(x_test)["f(x)"]
         sigma = self.posterior_covariance(x_test)["v(x)"]
-        r = self._crps_s(y_test, mean, sigma)
+        r = self._crps_s(y_test, mean, np.sqrt(sigma))
         return r
 
     def rmse(self, x_test, y_test):
@@ -1326,6 +1326,56 @@ def rmse(self, x_test, y_test):
         v2 = self.posterior_mean(x_test)["f(x)"].reshape(len(v1))
         return np.sqrt(np.sum((v1 - v2) ** 2) / len(v1))
 
+    def nlpd(self, x_test, y_test):
+        """
+        This function calculates the Negative log predictive density.
+        Note that in the multi-task setting the user should perform their
+        input point transformation beforehand.
+
+        Parameters
+        ----------
+        x_test : np.ndarray
+            A numpy array of shape (V x D), interpreted as  an array of input point positions.
+        y_test : np.ndarray
+            A numpy array of shape V. These are the y data to compare against
+
+        Return
+        ------
+        NLPD : float
+        """
+
+        mean = self.posterior_mean(x_test)["f(x)"]
+        sigma = np.sqrt(self.posterior_covariance(x_test)["v(x)"])
+
+        g = self.gaussian_1d(y_test, mean, sigma)
+        g[g == 0.] = 1e-16
+        g = np.log(g)
+        return -np.mean(g)
+
+    @staticmethod
+    def gaussian_1d(x, mu, sigma):
+        """
+        Evaluates a 1D Gaussian (Normal) distribution at a point x.
+
+        Parameters
+        ----------
+        x : np.ndarray
+            The points where you want to evaluate the Gaussian.
+        mu : np.ndarray
+            The mean of the Gaussian (default 0.0).
+        sigma : np.ndarray
+            The standard deviation of the Gaussians.
+
+        Return
+        ------
+        Evaluations of the Gaussian : np.ndarray
+        """
+        # Gaussian function formula
+        coefficient = 1.0 / (np.sqrt(2 * np.pi) * sigma)
+        exponent = -((x - mu) ** 2) / (2 * sigma ** 2)
+        return coefficient * np.exp(exponent)
+
+
     @staticmethod
     def make_2d_x_pred(bx, by, resx=100, resy=100):  # pragma: no cover
         """

tests/test_fvgp.py

@@ -116,6 +116,7 @@ def prior_mean(x,hps):
     res = my_gp1.prior._default_kernel(x_data,x_data,np.array([1.,1.,1.,1.,1.,1.]))
     my_gp1.crps(x_data[0:2] + 1., np.array([1.,2.]))
     my_gp1.rmse(x_data[0:2] + 1., np.array([1.,2.]))
+    my_gp1.nlpd(x_data[0:2] + 1., np.array([1.,2.]))
     my_gp1.make_2d_x_pred(np.array([1.,2.]),np.array([3.,4]))
     my_gp1.make_1d_x_pred(np.array([1.,2.]))
     my_gp1._get_default_hyperparameter_bounds()