数学建模习题7.10

`import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy.interpolate import interp1d, PchipInterpolator, CubicSpline
from scipy.optimize import curve_fit
from scipy.stats import norm

读取数据

file_path = '7.17.xISX' # 确保文件路径正确
data = pd.read_excel(file_path, header=None)
x_data = data.iloc[:, 0].values
y_data = data.iloc[:, 1].values

定义插值函数

x_interp = np.linspace(-2, 4.9, 400)
interp_functions = {
'Linear': interp1d(x_data, y_data, kind='linear', bounds_error=False, fill_value='extrapolate'),
'Cubic': interp1d(x_data, y_data, kind='cubic', bounds_error=False, fill_value='extrapolate'),
'PCHIP': PchipInterpolator(x_data, y_data, extrapolate=True),
'CubicSpline': CubicSpline(x_data, y_data, bc_type='natural', extrapolate=True)
}
y_interps = {name: func(x_interp) for name, func in interp_functions.items()}

绘制插值结果

plt.figure(figsize=(10, 6))
for name, y_interp in y_interps.items():
plt.plot(x_interp, y_interp, label=name)
plt.plot(x_data, y_data, 'o', label='Original Data')
plt.xlabel('x')
plt.ylabel('y')
plt.legend()
plt.title('Interpolation Methods Comparison')
plt.grid(True)
plt.show()

多项式拟合

degrees = range(1, 6)
coeffs = {}
residuals = {}
for degree in degrees:
coefficients, residual = np.polyfit(x_data, y_data, degree, cov=True)
residual_std = np.sqrt(residual[0, 0])
coeffs[degree] = coefficients
residuals[degree] = residual_std

for degree, coeffs_val in coeffs.items():
print(f"Degree {degree} Polynomial Coefficients: {coeffs_val}")
print(f"Residual Standard Deviation: {residuals[degree]:.4f}")

best_degree = min(residuals, key=residuals.get)
print(f"Best fitting polynomial degree: {best_degree}")
best_poly = np.poly1d(coeffs[best_degree])
y_poly_fit = best_poly(x_interp)

绘制多项式拟合结果

plt.figure(figsize=(10, 6))
plt.plot(x_interp, y_poly_fit, label=f"{best_degree}-degree Polynomial Fit")
plt.plot(x_data, y_data, 'o', label='Original Data')
plt.xlabel('x')
plt.ylabel('y')
plt.legend()
plt.title('Polynomial Fit')
plt.grid(True)
plt.show()

正态分布拟合

def normal_dist(x, mu, sigma):
return norm.pdf(x, mu, sigma)

params, params_covariance = curve_fit(normal_dist, x_data, y_data, p0=[np.mean(x_data), np.std(x_data)])
mu, sigma = params
x_normal_fit = np.linspace(min(x_data), max(x_data), 400)
y_normal_fit = normal_dist(x_normal_fit, mu, sigma)

绘制正态分布拟合结果

plt.figure(figsize=(10, 6))
plt.plot(x_normal_fit, y_normal_fit, label=f'Normal Distribution Fit\nmu={mu:.2f}, sigma={sigma:.2f}')
plt.plot(x_data, y_data, 'o', label='Original Data')
plt.xlabel('x')
plt.ylabel('y')
plt.legend()
plt.title('Normal Distribution Fit')
plt.grid(True)
plt.show()

print("学号后四位：3005")`