Practical Applications & Tutorials#
This section provides detailed examples of using torch-secorder in various scenarios. For basic usage of foundational tools like Hessian-Vector Products (HVPs), Jacobian-Vector Products (JVPs), and basic curvature statistics (Hessian diagonal/trace), please refer to their respective API documentation in the Core and Approximations modules.
Analysis Examples#
Loss Landscape Visualization#
This example demonstrates how to compute and visualize the loss landscape of a neural network model. The loss landscape provides insights into the geometry of the optimization space and can help understand:
The smoothness and convexity of the loss function
The presence of local minima and saddle points
The difficulty of optimization in different regions
The example shows two types of visualizations:
1D Loss Surface: A slice of the loss landscape along a random direction, showing how the loss changes as we move away from the current parameter point.
2D Loss Surface: A contour plot showing the loss landscape in a 2D plane defined by two random directions, providing a more comprehensive view of the optimization space.
The visualization generates two files: - 1d_loss_surface.png: A line plot showing the loss along a single direction - 2d_loss_surface.png: A 3D surface plot showing the loss in a 2D plane
1"""Example demonstrating loss landscape visualization.
2
3This script shows how to compute and visualize 1D slices and 2D contours
4of the loss landscape for a simple model, using random directions.
5
6The example demonstrates:
71. Computing 1D loss surface along a random direction
82. Computing 2D loss surface using two random directions
93. Visualizing both surfaces using matplotlib
10
11Requirements:
12 - torch
13 - matplotlib
14 - numpy
15"""
16
17import matplotlib.pyplot as plt
18import torch
19import torch.nn as nn
20import torch.nn.functional as F
21from matplotlib import cm
22
23from torch_secorder.analysis.landscape import (
24 compute_loss_surface_1d,
25 compute_loss_surface_2d,
26 create_random_direction,
27)
28
29
30# 1. Define a Simple Model
31class SimpleModel(nn.Module):
32 def __init__(self):
33 super().__init__()
34 self.linear = nn.Linear(2, 1)
35
36 def forward(self, x):
37 return self.linear(x)
38
39
40# 2. Generate Synthetic Data
41torch.manual_seed(42)
42X_train = torch.randn(20, 2)
43y_train = X_train @ torch.tensor([[2.0], [3.0]]) + 1.0 + torch.randn(20, 1) * 0.5
44
45
46# 3. Instantiate Model and Loss Function
47model = SimpleModel()
48loss_fn = F.mse_loss
49
50
51print("--- Loss Landscape Visualization Example ---")
52
53# 4. Compute 1D Loss Surface
54print("\nComputing 1D loss surface...")
55direction_1d = create_random_direction(model)
56alphas_1d, losses_1d = compute_loss_surface_1d(
57 model,
58 loss_fn,
59 X_train,
60 y_train,
61 direction_1d,
62 alpha_range=(-2.0, 2.0),
63 num_points=50,
64)
65
66print("1D Loss Surface (first 5 points):\nAlphas: ", alphas_1d[:5].tolist())
67print("Losses: ", losses_1d[:5].tolist())
68
69# Plot 1D Loss Surface
70plt.figure(figsize=(8, 6))
71plt.plot(alphas_1d.numpy(), losses_1d.numpy(), "b-", linewidth=2)
72plt.xlabel("Alpha (Direction Scale)")
73plt.ylabel("Loss")
74plt.title("1D Loss Surface")
75plt.grid(True)
76plt.savefig("1d_loss_surface.png")
77plt.close()
78
79
80# 5. Compute 2D Loss Surface
81print("\nComputing 2D loss surface...")
82direction1_2d = create_random_direction(model)
83direction2_2d = create_random_direction(model)
84
85# Ensure directions are not collinear (optional, but good for meaningful 2D surface)
86# A simple way to get somewhat orthogonal directions: re-randomize if dot product is too high
87# This is a heuristic, proper orthogonalization methods might be preferred for robustness
88if (
89 torch.dot(
90 torch.cat([d.flatten() for d in direction1_2d]),
91 torch.cat([d.flatten() for d in direction2_2d]),
92 ).abs()
93 > 0.5
94):
95 print("Adjusting second random direction for better orthogonality...")
96 direction2_2d = create_random_direction(model)
97
98alphas_2d, betas_2d, losses_2d = compute_loss_surface_2d(
99 model, loss_fn, X_train, y_train, direction1_2d, direction2_2d, num_points=25
100)
101
102print("2D Loss Surface (top-left 3x3 values):\n", losses_2d[:3, :3].tolist())
103
104# Plot 2D Loss Surface
105fig = plt.figure(figsize=(12, 10))
106ax = fig.add_subplot(111, projection="3d")
107A, B = torch.meshgrid(alphas_2d, betas_2d, indexing="ij")
108surf = ax.plot_surface(
109 A.numpy(),
110 B.numpy(),
111 losses_2d.numpy(),
112 cmap=cm.viridis,
113 edgecolor="none",
114 alpha=0.8,
115)
116ax.set_xlabel("Alpha (Direction 1)")
117ax.set_ylabel("Beta (Direction 2)")
118ax.set_zlabel("Loss")
119ax.set_title("2D Loss Surface")
120fig.colorbar(surf, shrink=0.5, aspect=5)
121plt.savefig("2d_loss_surface.png")
122plt.close()
123
124print(
125 "\nLoss landscape visualization complete. Check '1d_loss_surface.png' and '2d_loss_surface.png' for the plots."
126)
Per-Layer Curvature Analysis#
This example demonstrates how to extract per-layer Fisher diagonal elements, which form the basis for Kronecker-Factored Approximate Curvature (K-FAC) methods.
1"""Example demonstrating per-layer curvature extraction and its relation to K-FAC.
2
3This script shows how to extract per-layer Fisher diagonal elements using per_layer_fisher_diagonal.
4These per-layer blocks are the foundation for Kronecker-Factored Approximate Curvature (K-FAC).
5"""
6
7import torch
8import torch.nn as nn
9import torch.nn.functional as F
10
11from torch_secorder.core.per_layer_curvature import per_layer_fisher_diagonal
12
13
14# Define a simple model with multiple layers
15class SimpleModel(nn.Module):
16 def __init__(self):
17 super().__init__()
18 self.linear1 = nn.Linear(10, 5)
19 self.linear2 = nn.Linear(5, 3)
20
21 def forward(self, x):
22 x = F.relu(self.linear1(x))
23 return self.linear2(x)
24
25
26# Generate synthetic data
27torch.manual_seed(42)
28X = torch.randn(8, 10)
29y = torch.randint(0, 3, (8,)) # 3 classes
30
31# Instantiate the model
32model = SimpleModel()
33
34# Compute per-layer Fisher diagonal
35# This extracts block-diagonal approximations of the Fisher Information Matrix
36# These blocks are the basis for K-FAC approximations
37layer_fisher_diagonals = per_layer_fisher_diagonal(model, F.cross_entropy, X, y)
38
39print("Per-Layer Fisher Diagonal (K-FAC basis):")
40for layer_name, diag in layer_fisher_diagonals.items():
41 print(
42 f"Layer: {layer_name}, Shape: {diag.shape}, First few values: {diag.flatten()[:5].tolist()}"
43 )
44
45print(
46 "\nK-FAC would further decompose these blocks into Kronecker products of smaller matrices."
47)
Optimization Examples#
Natural Gradient Descent#
This example demonstrates a simplified implementation of Natural Gradient Descent, leveraging the Fisher diagonal approximation from torch-secorder.
1"""Example demonstrating Natural Gradient Descent using Fisher diagonal approximation."""
2
3import torch
4import torch.nn as nn
5import torch.optim as optim
6
7from torch_secorder.approximations.fisher_diagonal import empirical_fisher_diagonal
8
9
10# 1. Define a Simple Model
11class SimpleModel(nn.Module):
12 def __init__(self):
13 super().__init__()
14 self.linear = nn.Linear(2, 1)
15
16 def forward(self, x):
17 return self.linear(x)
18
19
20# 2. Generate Synthetic Data
21torch.manual_seed(42)
22X = torch.randn(100, 2)
23y = X @ torch.tensor([[1.0], [2.0]]) + 0.5 + torch.randn(100, 1) * 0.1
24
25# 3. Instantiate Model, Loss, and Optimizer
26model = SimpleModel()
27loss_fn = nn.MSELoss()
28
29# In a real NGD optimizer, this would be integrated directly.
30# Here, we demonstrate the manual step using Fisher diagonal.
31
32learning_rate = 0.01
33num_epochs = 100
34
35print("\n--- Training with Natural Gradient Descent (Approximation) ---")
36for epoch in range(num_epochs):
37 model.train()
38 optimizer = optim.SGD(
39 model.parameters(), lr=learning_rate
40 ) # We'll manually adjust gradients
41 optimizer.zero_grad()
42
43 # Forward pass
44 outputs = model(X)
45 loss = loss_fn(outputs, y)
46
47 # Backward pass to compute gradients
48 loss.backward(
49 create_graph=True
50 ) # create_graph=True is needed for higher-order derivatives like Fisher
51
52 # Compute Empirical Fisher Diagonal
53 # For NGD, we need F^{-1}g. Here, we approximate F as diagonal.
54 fisher_diagonal_list = empirical_fisher_diagonal(model, loss_fn, X, y)
55
56 # Manual Natural Gradient Update (approximated with diagonal Fisher)
57 with torch.no_grad():
58 for i, param in enumerate(model.parameters()):
59 if param.grad is not None:
60 # Natural gradient: F^{-1} * gradient
61 # If F is diagonal, F_ii = fisher_diagonal_list[i]_ii
62 # So F^{-1}_ii = 1 / fisher_diagonal_list[i]_ii
63 # Natural gradient update for each parameter:
64 # param_update = learning_rate * (param.grad / fisher_diagonal_list[i])
65 # param.data -= param_update
66 # A simpler diagonal update often used in practice for NGD is param.grad / (F_diag + epsilon)
67
68 # Let's use a simpler diagonal update: param.grad / (fisher_diag + epsilon)
69 epsilon = 1e-6 # For numerical stability
70 param.grad.data.div_(
71 fisher_diagonal_list[i].clamp(min=epsilon)
72 ) # Element-wise division
73
74 optimizer.step() # Apply the modified gradients
75
76 if (epoch + 1) % 10 == 0:
77 print(f"Epoch [{epoch + 1}/{num_epochs}], Loss: {loss.item():.4f}")
78
79print("\nTraining complete.")
80print(f"Final Loss: {loss.item():.4f}")
81print("Model parameters after NGD:")
82for name, param in model.named_parameters():
83 print(f"{name}: {param.data.squeeze()}")
Bayesian Neural Networks Inference#
This example illustrates how Fisher Information can be used in the context of Bayesian Neural Networks to approximate posterior variances of parameters.
1"""Example demonstrating the use of Fisher Information for Bayesian Neural Networks (BNNs)."""
2
3import torch
4import torch.nn as nn
5import torch.nn.functional as F
6
7from torch_secorder.approximations import (
8 empirical_fisher_diagonal,
9 generalized_fisher_diagonal,
10)
11
12
13# 1. Define a Simple Bayesian Neural Network (for demonstration purposes)
14# In a true BNN, weights would be distributions, but here we'll use point estimates
15# and use Fisher to infer uncertainty or for approximate posterior.
16class SimpleBNN(nn.Module):
17 def __init__(self):
18 super().__init__()
19 self.linear1 = nn.Linear(10, 5)
20 self.linear2 = nn.Linear(5, 1)
21
22 def forward(self, x):
23 x = F.relu(self.linear1(x))
24 return self.linear2(x)
25
26
27# 2. Generate Synthetic Data
28torch.manual_seed(42)
29X_train = torch.randn(100, 10)
30y_train = X_train @ torch.randn(10, 1) + 0.5 + torch.randn(100, 1) * 0.1
31
32X_test = torch.randn(20, 10)
33
34# 3. Instantiate Model and Loss Function
35model = SimpleBNN()
36
37# For BNNs, often working with negative log-likelihood (NLL) directly
38# For regression, this might be Gaussian NLL. For classification, CrossEntropy is NLL.
39# Here, we'll use MSE for simplicity and demonstrate Fisher diagonal for it.
40# Note: Generalized Fisher is typically for likelihoods (e.g., NLL in classification).
41# For MSE in regression, Empirical Fisher is often more directly applicable.
42
43# Let's simulate a regression task with a Gaussian likelihood assumption
44# where the loss is proportional to MSE (which is -log-likelihood for Gaussian with fixed variance)
45# We'll adapt `generalized_fisher_diagonal` for a simplified regression NLL scenario.
46
47print("\n--- Estimating Uncertainty using Fisher Information (BNN Context) ---")
48
49# Train the model (simplified training loop)
50optimizer = torch.optim.Adam(model.parameters(), lr=0.01)
51for epoch in range(50):
52 optimizer.zero_grad()
53 outputs = model(X_train)
54 loss = F.mse_loss(outputs, y_train)
55 loss.backward()
56 optimizer.step()
57 if (epoch + 1) % 10 == 0:
58 print(f"Epoch [{epoch + 1}/50], Loss: {loss.item():.4f}")
59
60# 4. Use Fisher Information for Uncertainty Estimation (Approximation)
61# For a BNN, the inverse Fisher (or Hessian) often approximates the posterior covariance.
62# The diagonal of F^{-1} gives approximate posterior variances of parameters.
63
64print(
65 "\n--- Computing Fisher Diagonal for Approximate Posterior Variance (Regression: MSE) ---"
66)
67
68print(
69 "\nUsing Empirical Fisher Diagonal for parameter uncertainty approximation (MSE regression)..."
70)
71try:
72 # Using empirical_fisher_diagonal for a direct demonstration applicable to MSE
73 fisher_diagonal_params = empirical_fisher_diagonal(
74 model, F.mse_loss, X_train, y_train
75 )
76
77 print("\nApproximate Posterior Variances (from inverse Fisher diagonal):")
78 with torch.no_grad():
79 for i, param in enumerate(model.parameters()):
80 # Inverse of Fisher diagonal approximates variance of parameters
81 # Ensure no division by zero if Fisher diagonal element is too small
82 variance_approx = 1.0 / (
83 fisher_diagonal_params[i] + 1e-8
84 ) # Add epsilon for stability
85 print(
86 f"Parameter {i} variance approximation shape: {variance_approx.shape}"
87 )
88 print(f" First few values: {variance_approx.flatten()[:5].tolist()}")
89
90except NotImplementedError as e:
91 print(
92 f"Error: {e}. Generalized Fisher diagonal for regression NLL might not be fully supported yet in `generalized_fisher_diagonal` for direct MSE loss."
93 )
94except Exception as e:
95 print(f"An unexpected error occurred: {e}")
96
97print("\nBNN Fisher Information example complete.")
98
99# Add example for generalized_fisher_diagonal (classification use case)
100print("\n--- Using Generalized Fisher Diagonal (Classification Example) ---")
101
102
103# Create a simple classification model and data
104class SimpleClassifier(nn.Module):
105 def __init__(self):
106 super().__init__()
107 self.linear = nn.Linear(10, 3)
108
109 def forward(self, x):
110 return self.linear(x)
111
112
113clf_model = SimpleClassifier()
114clf_inputs = torch.randn(8, 10)
115clf_targets = torch.randint(0, 3, (8,)) # 3 classes
116clf_outputs = clf_model(clf_inputs)
117
118# Compute the Generalized Fisher Diagonal (NLL loss type)
119gfim_diagonal = generalized_fisher_diagonal(
120 clf_model, clf_outputs, clf_targets, loss_type="nll"
121)
122print("Generalized Fisher Diagonal (classification):")
123for i, diag in enumerate(gfim_diagonal):
124 print(f" Param {i} shape: {diag.shape}, first few: {diag.flatten()[:5].tolist()}")
Gauss-Newton Optimization#
This example demonstrates a simplified optimization loop using the diagonal Gauss-Newton Matrix approximation.
1"""Example demonstrating a simplified Gauss-Newton optimization step using GNM diagonal."""
2
3import torch
4import torch.nn as nn
5import torch.optim as optim
6
7from torch_secorder.approximations.gauss_newton import gauss_newton_matrix_approximation
8
9
10# 1. Define a Simple Model (e.g., for non-linear least squares)
11class SimpleNonLinearModel(nn.Module):
12 def __init__(self):
13 super().__init__()
14 # Let's make it slightly non-linear
15 self.linear1 = nn.Linear(2, 5)
16 self.relu = nn.ReLU()
17 self.linear2 = nn.Linear(5, 1)
18
19 def forward(self, x):
20 return self.linear2(self.relu(self.linear1(x)))
21
22
23# 2. Generate Synthetic Data
24torch.manual_seed(42)
25X = torch.randn(100, 2) * 5 # Input features
26y = (X[:, 0] ** 2 + X[:, 1] * 3 + 2 + torch.randn(100) * 0.5).unsqueeze(
27 1
28) # Non-linear target
29
30# 3. Instantiate Model, Loss
31model = SimpleNonLinearModel()
32loss_fn = nn.MSELoss() # Gauss-Newton is for least squares
33
34learning_rate = 0.01
35num_iterations = 100
36
37print("\n--- Training with Approximate Gauss-Newton (Diagonal GNM) ---")
38for iteration in range(num_iterations):
39 model.train()
40 optimizer = optim.SGD(
41 model.parameters(), lr=learning_rate
42 ) # We'll manually adjust gradients
43 optimizer.zero_grad()
44
45 # Forward pass
46 outputs = model(X)
47 loss = loss_fn(outputs, y)
48
49 # Backward pass for gradients
50 loss.backward(
51 create_graph=True
52 ) # create_graph=True if we want higher-order, not strictly needed for GNM diag
53
54 # Compute Diagonal of Gauss-Newton Matrix
55 # The inverse of the GNM is used in the update rule: delta_theta = - G^{-1} * gradient
56 # Here, we use diagonal GNM, so G^{-1}_ii = 1 / G_ii
57 gnm_diagonal_list = list(gauss_newton_matrix_approximation(model, loss_fn, X, y))
58
59 # Manual Gauss-Newton Update (approximated with diagonal GNM)
60 with torch.no_grad():
61 for i, param in enumerate(model.parameters()):
62 if param.grad is not None:
63 # Update: param.data -= learning_rate * (param.grad / (gnm_diag + epsilon))
64 epsilon = 1e-6 # For numerical stability
65 # Ensure gnm_diagonal_list[i] has the same shape as param.grad
66 param.grad.data.div_(
67 gnm_diagonal_list[i].clamp(min=epsilon)
68 ) # Element-wise division
69
70 optimizer.step() # Apply the modified gradients
71
72 if (iteration + 1) % 10 == 0:
73 print(f"Iteration [{iteration + 1}/{num_iterations}], Loss: {loss.item():.4f}")
74
75print("\nTraining complete.")
76print(f"Final Loss: {loss.item():.4f}")
77print("Model parameters after Gauss-Newton approximation:")
78for name, param in model.named_parameters():
79 print(f"{name}: {param.data.squeeze()}")
General Examples#
Quick Start#
A minimal example showing how to get started with torch-secorder.
1"""Quick start example for torch-secorder.
2
3This example demonstrates basic usage of torch-secorder's core functionality:
41. Computing Hessian-Vector Products (HVPs)
52. Computing Jacobian-Vector Products (JVPs)
63. Computing Hessian diagonal
7"""
8
9import torch
10import torch.nn as nn
11import torch.nn.functional as F
12
13from torch_secorder.core import gvp, hessian_diagonal, hvp
14
15
16def main():
17 # 1. Create a simple model
18 model = nn.Sequential(nn.Linear(10, 5), nn.ReLU(), nn.Linear(5, 1))
19
20 # 2. Generate some random data
21 x = torch.randn(32, 10)
22 y = torch.randn(32, 1)
23
24 # 3. Compute loss
25 def get_loss():
26 return F.mse_loss(model(x), y)
27
28 # 4. Create a random direction vector
29 v = [torch.randn_like(p) for p in model.parameters()]
30
31 # 5. Compute HVP
32 hvp_result = hvp.model_hvp(model, F.mse_loss, x, y, v)
33 print("HVP computed successfully!")
34
35 # 6. Compute JVP
36 jvp_result = gvp.model_jvp(model, x, v)
37 print("JVP computed successfully!")
38
39 # 7. Compute Hessian diagonal
40 hessian_diag = hessian_diagonal(get_loss, list(model.parameters()))
41 print("Hessian diagonal computed successfully!")
42
43 # 8. Print shapes for verification
44 print("\nShapes of results:")
45 print(f"HVP: {[h.shape for h in hvp_result]}")
46 print(f"JVP: {[j.shape for j in jvp_result]}")
47 print(f"Hessian diagonal: {[h.shape for h in hessian_diag]}")
48
49
50if __name__ == "__main__":
51 main()
Custom Model Integration#
Demonstrates how to integrate torch-secorder with custom PyTorch models.
1"""Example demonstrating integration with custom PyTorch models.
2
3This example shows how to use torch-secorder with custom PyTorch models,
4including models with:
51. Custom forward passes
62. Multiple outputs
73. Complex loss functions
8"""
9
10import torch
11import torch.nn as nn
12import torch.nn.functional as F
13
14from torch_secorder.core import gvp, hessian_diagonal, hvp
15
16
17class CustomModel(nn.Module):
18 """A custom model with multiple outputs and a complex forward pass."""
19
20 def __init__(self):
21 super().__init__()
22 self.encoder = nn.Sequential(nn.Linear(10, 20), nn.ReLU(), nn.Linear(20, 10))
23 self.classifier = nn.Linear(10, 2)
24 self.regressor = nn.Linear(10, 1)
25
26 def forward(self, x):
27 features = self.encoder(x)
28 return {
29 "logits": self.classifier(features),
30 "predictions": self.regressor(features),
31 }
32
33
34def custom_loss(outputs, targets):
35 """A custom loss function combining classification and regression losses."""
36 classification_loss = F.cross_entropy(outputs["logits"], targets["classes"])
37 regression_loss = F.mse_loss(outputs["predictions"], targets["values"])
38 return classification_loss + 0.1 * regression_loss
39
40
41def main():
42 # 1. Create the custom model
43 model = CustomModel()
44
45 # 2. Generate random data
46 x = torch.randn(32, 10)
47 y = {"classes": torch.randint(0, 2, (32,)), "values": torch.randn(32, 1)}
48
49 # 3. Compute loss
50 def get_loss():
51 return custom_loss(model(x), y)
52
53 # 4. Create a random direction vector
54 v = [torch.randn_like(p) for p in model.parameters()]
55
56 # 5. Compute HVP with custom model and loss
57 hvp_result = hvp.model_hvp(model, custom_loss, x, y, v)
58 print("HVP computed successfully with custom model!")
59
60 # 6. Compute JVP for each output separately
61 # For models with multiple outputs, we need to compute JVP for each output tensor
62 # We'll create wrapper models for each output
63
64 class LogitsModel(nn.Module):
65 def __init__(self, base_model):
66 super().__init__()
67 self.base_model = base_model
68
69 def forward(self, x):
70 return self.base_model(x)["logits"]
71
72 class PredictionsModel(nn.Module):
73 def __init__(self, base_model):
74 super().__init__()
75 self.base_model = base_model
76
77 def forward(self, x):
78 return self.base_model(x)["predictions"]
79
80 # Create wrapper models
81 logits_model = LogitsModel(model)
82 predictions_model = PredictionsModel(model)
83
84 # Compute JVPs using the wrapper models
85 jvp_logits = gvp.model_jvp(logits_model, x, v)
86 print("JVP computed successfully for logits!")
87
88 jvp_predictions = gvp.model_jvp(predictions_model, x, v)
89 print("JVP computed successfully for predictions!")
90
91 # 7. Compute Hessian diagonal
92 hessian_diag = hessian_diagonal(get_loss, list(model.parameters()))
93 print("Hessian diagonal computed successfully with custom model!")
94
95 # 8. Print shapes for verification
96 print("\nShapes of results:")
97 print(f"HVP: {[h.shape for h in hvp_result]}")
98 print(f"JVP (logits): {[j.shape for j in jvp_logits]}")
99 print(f"JVP (predictions): {[j.shape for j in jvp_predictions]}")
100 print(f"Hessian diagonal: {[h.shape for h in hessian_diag]}")
101
102
103if __name__ == "__main__":
104 main()