申し訳ありません。スクリプトの内容を直接お伝えできていませんでした。以下に `generate_acquisition_plots.py` のコード内容を提供します。
```python
import numpy as np
import matplotlib.pyplot as plt
from skopt import gp_minimize
from skopt.space import Real
from skopt.learning import GaussianProcessRegressor
from skopt.utils import create_acquisition_function
# Optimization history
X_history = [[5.0], [3.07], [3.15]]
y_history = [48.1852, 28.773, 29.5509]
P_target = 29.589892
# Objective function (dummy, as we already have results)
def objective(U_wall):
# This function is not called in this script for actual evaluation
# It is merely a placeholder for the skopt interface
pass
# Define the search space
space = [Real(0.1, 10.0, name='U_wall')]
# Function to plot acquisition function
def plot_acquisition(X, y, gp_model, iteration):
fig, ax1 = plt.subplots(figsize=(10, 6))
# Plot observed points
ax1.scatter(np.array(X).flatten(), y, c='red', label='Observed points')
# Create a dense grid for plotting
x_plot = np.linspace(0.1, 10.0, 1000).reshape(-1, 1)
# Predict mean and std for the dense grid
mu, sigma = gp_model.predict(x_plot, return_std=True)
ax1.plot(x_plot, mu, label='Mean prediction')
ax1.fill_between(x_plot.flatten(), mu - 1.96 * sigma, mu + 1.96 * sigma, alpha=0.2, label='95% confidence interval')
ax1.axhline(P_target, color='green', linestyle='--', label='Target Pressure')
ax1.set_xlabel('U_wall (m/s)')
ax1.set_ylabel('Pressure Difference')
ax1.tick_params(axis='y', labelcolor='blue')
ax1.set_title(f'Iteration {iteration} - Surrogate Model and Acquisition Function')
ax1.legend(loc='upper left')
# Create a secondary y-axis for the acquisition function
ax2 = ax1.twinx()
# Use Expected Improvement (EI) for acquisition function
# Note: skopt's acquisition functions are usually designed to minimize, so we might need to adjust or create a custom one if directly targeting P_target
# For this example, we'll use the default acquisition function (expected_improvement) with the adjusted objective to minimize the absolute difference from target
# A custom acquisition function can be built, but for simplicity, let's demonstrate the built-in EI
# We are minimizing the absolute difference from the target
y_adjusted = np.abs(np.array(y) - P_target)
# Re-train GP model for acquisition function calculation based on adjusted y
gp_model_adjusted = GaussianProcessRegressor(alpha=1e-10, normalize_y=True)
gp_model_adjusted.fit(X, y_adjusted)
acq_func = create_acquisition_function(gp_model_adjusted, y_adjusted, "EI")
acquisition_values = acq_func(x_plot)
ax2.plot(x_plot, acquisition_values, color='purple', linestyle='-', label='Acquisition Function (EI)')
ax2.set_ylabel('Acquisition Function Value', color='purple')
ax2.tick_params(axis='y', labelcolor='purple')
ax2.legend(loc='upper right')
# Find the next recommended point (maximizer of the acquisition function)
next_x_idx = np.argmax(acquisition_values)
next_x = x_plot[next_x_idx][0]
ax2.axvline(next_x, color='blue', linestyle='--', label=f'Next U_wall: {next_x:.2f} m/s')
ax2.legend(loc='upper right')
plt.grid(True)
plt.tight_layout()
plt.savefig(f'acquisition_plot_iteration_{iteration}.png')
plt.close(fig)
# Simulate the optimization steps and plot
for i in range(1, len(X_history) + 1):
if i > 1: # Only plot after at least one observation for meaningful GP
X_current = X_history[:i]
y_current = y_history[:i]
# Initialize and fit Gaussian Process Regressor
gp_model = GaussianProcessRegressor(alpha=1e-10, normalize_y=True)
gp_model.fit(X_current, y_current)
# Plot for the current iteration
plot_acquisition(X_current, y_current, gp_model, i)
print("Acquisition function plots generated: acquisition_plot_iteration_2.png, acquisition_plot_iteration_3.png")
```
