最小二乘法实际应用：大气二氧化碳浓度数据拟合

最小二乘法实际应用：大气二氧化碳浓度数据拟合

引用

CSDN

1.

https://m.blog.csdn.net/qq_50930131/article/details/144274004

最小二乘法是一种数学优化技术，它通过最小化误差的平方和寻找数据的最佳函数匹配。这种方法广泛应用于回归分析中，特别是在线性回归模型中。本文将通过一个实际案例，演示如何使用最小二乘法拟合大气二氧化碳浓度数据。

使用最小二乘法拟合大气二氧化碳浓度数据。数据保存在monthly_co2.xls文件中（只截取部分）。

Python需要安装的库

1. xlrd
2. numpy
3. pandas
4. matplotlib

绘制图像代码（绘制整体数据趋势图）

# -*- coding: utf-8 -*- 
"""
@File    :   绘制趋势图.py
@Time    :   2024/11/27 23:52:44
@Version :   
@Desc    :   
"""
"""
@Python version         : 3.8.7
@matplotlib version     : 3.5.1
@pandas version         : 1.4.2
@numpy version          : 1.22.3
"""
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from matplotlib.ticker import MultipleLocator                   #从ticker中导入MultipleLocator类
df = pd.read_excel(io="./monthly_co2.xls", header=2)            # 数据在第一个sheet中
# 删除索引为0的行
df.drop(index=0, inplace=True)                                  # 删除第三行数据    直接在原始dataframe上修改   
df.replace(-99.99, np.nan, inplace=True)                        # 处理缺失值
df.dropna(inplace=True)                                         # 删除缺失值的行
# df['Year_Month'] = df['Yr'] + (df['Mn'] - 1) / 12             # 将年份和月份转换为小数形式
# [:, 0]    表示选择所有行
df['Year_Month'] = df.iloc[:, 0] + (df.iloc[:, 1] - 1) / 12     # 使用列index索引而不是名称索引
x_vec = df['Year_Month'].values                                 # x
y_vec = df.iloc[:, 2].values                                    # y
# matplotlib字体设置
plt.rcParams['font.family'] = "Times New Roman"                 # 设置全局字体
# marker='o'    散点图绘制为圆形
# edgecolor='r' 设置圆圈的边缘颜色为红色
custom_color = '#f0a1a8'
plt.scatter(x_vec, y_vec, marker='o', edgecolors=custom_color, facecolors='none', label='Actual Data')                  # 原始数据
plt.xlabel("Year_Month", fontsize=16)
plt.ylabel("CO2 Concentration [ppm]", fontsize=16)
plt.title("Trend Chart of CO2 Concentration Over Time", fontsize=20)         # 二氧化碳浓度随时间变化趋势图
plt.legend(loc='upper left')
# 美化图表
plt.rcParams['axes.facecolor'] = 'lightgray'
plt.rcParams['legend.fontsize'] = 10                            # 设置图例字体大小
# 刻度值设置
ax = plt.gca()                                                  # 获取当前坐标轴
ax.xaxis.set_major_locator(MultipleLocator(10))                 # 设置x轴的主刻度间隔
ax.xaxis.set_minor_locator(MultipleLocator(5))                  # 设置x轴的次刻度间隔
ax.yaxis.set_major_locator(MultipleLocator(20))                 # 设置y轴的主刻度间隔
ax.yaxis.set_minor_locator(MultipleLocator(10))                 # 设置y轴的次刻度间隔
plt.grid(True, color='gray', linestyle='--', linewidth=0.5)     # 添加自定义样式网格线
plt.show()

二次函数拟合

# -*- coding: utf-8 -*- 
"""
@File    :   二次函数拟合结果.py
@Time    :   2024/11/28 11:19:38
@Version :   
@Desc    :   
"""
"""
@Python version         : 3.8.7
@matplotlib version     : 3.5.1
@pandas version         : 1.4.2
@numpy version          : 1.22.3
"""
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from matplotlib.ticker import MultipleLocator                   #从ticker中导入MultipleLocator类
df = pd.read_excel(io="./monthly_co2.xls", header=2)            # 数据在第一个sheet中
![](https://wy-static.wenxiaobai.com/chat-rag-image/16219254573352230935)
# 删除索引为0的行
df.drop(index=0, inplace=True)                                  # 删除第三行数据    直接在原始dataframe上修改   
df.replace(-99.99, np.nan, inplace=True)                        # 处理缺失值
df.dropna(inplace=True)                                         # 删除缺失值的行
# df['Year_Month'] = df['Yr'] + (df['Mn'] - 1) / 12             # 将年份和月份转换为小数形式
# [:, 0]    表示选择所有行
df['Year_Month'] = df.iloc[:, 0] + (df.iloc[:, 1] - 1) / 12     # 使用列index索引而不是名称索引
x_vec = df['Year_Month'].values                                 # x
y_vec = df.iloc[:, 2].values                                    # y
# 这里添加代码
M = df.shape[0]                                                 # 数据点个数
# 方程组的系数矩阵
Phi = np.zeros((M, 3))                                          # 创建空矩阵
for i in range(M):
    Phi[i][0] = 1
    Phi[i][1] = x_vec[i]
    Phi[i][2] = x_vec[i]**2
# 正规方程组
A = Phi.T @ Phi
b = Phi.T @ y_vec
x = np.linalg.solve(A, b)
print(x)
a_0 = x[0]
a_1 = x[1]
a_2 = x[2]
# f_x = a_2x^2 + a_1x + a_0
varphi_x = lambda x: a_2 * x**2 + a_1 * x + a_0                        # 拟合函数
# 可视化图像
y_vec_model = varphi_x(x_vec)
# 均方根误差
RMSE = np.sqrt(np.sum((y_vec - y_vec_model)**2) / M)
print("均方根误差: ", RMSE)
# matplotlib字体设置
plt.rcParams['font.family'] = "Times New Roman"                 # 设置全局字体
plt.plot(x_vec, y_vec_model-279, label='Fitted Line', color='red')
# marker='o'    散点图绘制为圆形
# edgecolor='r' 设置圆圈的边缘颜色为红色
custom_color = '#f0a1a8'
plt.scatter(x_vec, y_vec-279, marker='o', edgecolors=custom_color, facecolors='none', label='Actual Data')                  # 原始数据
plt.xlabel("Year_Month", fontsize=16)
plt.ylabel("CO2 Concentration [ppm]", fontsize=16)
plt.title("Trend Chart of CO2 Concentration Over Time", fontsize=20)         # 二氧化碳浓度随时间变化趋势图
plt.legend(loc='upper left')
# 美化图表
plt.rcParams['axes.facecolor'] = 'lightgray'
plt.rcParams['legend.fontsize'] = 10                            # 设置图例字体大小
# 刻度值设置
ax = plt.gca()                                                  # 获取当前坐标轴
ax.xaxis.set_major_locator(MultipleLocator(10))                 # 设置x轴的主刻度间隔
ax.xaxis.set_minor_locator(MultipleLocator(5))                  # 设置x轴的次刻度间隔
ax.yaxis.set_major_locator(MultipleLocator(20))                 # 设置y轴的主刻度间隔
ax.yaxis.set_minor_locator(MultipleLocator(10))                 # 设置y轴的次刻度间隔
plt.grid(True, color='gray', linestyle='--', linewidth=0.5)     # 添加自定义样式网格线
# 添加RMSE文本信息
# (0.02, 0.85) 以轴的百分比来表示 分别对应x轴和y轴的为转移
plt.text(0.09, 0.95, f"RMSE: {RMSE:.2f}", transform=plt.gca().transAxes, fontsize=18)
plt.show()

程序运行结果：