基于Seaborn在折线拐点图中绘制两条最佳拟合线并求交点

核心思路 #

要实现无需人工选点的两段LOBF拟合，核心是自动识别数据拐点，将数据分割为陡降段与平台段后分别做线性回归。推荐用分段线性回归最小化残差法，比斜率阈值法更稳健。

具体实现步骤 #

1. 自动识别拐点与分段拟合 #

用scipy.optimize的曲线拟合功能，定义分段线性函数求解拐点及两段斜率：

import numpy as np
from scipy.optimize import curve_fit

# 定义分段线性函数
def piecewise_linear(x, x0, y0, k1, k2):
    return np.piecewise(x, [x < x0, x >= x0], 
                        [lambda x: k1*(x-x0)+y0, lambda x: k2*(x-x0)+y0])

# 代入数据拟合，popt返回拐点(x0,y0)、陡降斜率k1、平台斜率k2
x = np.array(your_x_data)
y = np.array(your_y_data)
# 可加初始参数p0提升拟合稳定性，比如p0=[x.mean(), y.mean(), -1, 0]
popt, _ = curve_fit(piecewise_linear, x, y)
x0, y0, k1, k2 = popt

2. 生成两段拟合线数据 #

根据拟合参数生成两段线的坐标点：

# 陡降段（x <= 拐点x0）
x_segment1 = x[x <= x0]
y_segment1 = k1*(x_segment1 - x0) + y0

# 平台段（x >= 拐点x0）
x_segment2 = x[x >= x0]
y_segment2 = k2*(x_segment2 - x0) + y0

3. Seaborn可视化与跨年度趋势展示 #

用sns.scatterplot展示原始数据，叠加拟合线，再通过FacetGrid实现跨年度对比：

import seaborn as sns
import matplotlib.pyplot as plt

# 单图展示分段拟合结果
sns.scatterplot(x='x_col', y='y_col', hue='year', data=your_df)
plt.plot(x_segment1, y_segment1, 'r-', label='陡降段LOBF')
plt.plot(x_segment2, y_segment2, 'b-', label='平台段LOBF')
plt.scatter(x0, y0, c='black', s=100, label='交点/拐点')
plt.legend()

# 跨年度FacetGrid展示
g = sns.FacetGrid(your_df, col='year', col_wrap=2)
g.map_dataframe(sns.scatterplot, x='x_col', y='y_col')

# 为每个年度子图添加分段拟合线
for ax, year in zip(g.axes, your_df['year'].unique()):
    year_data = your_df[your_df['year'] == year]
    x_year = year_data['x_col'].values
    y_year = year_data['y_col'].values
    popt_year, _ = curve_fit(piecewise_linear, x_year, y_year)
    x0_year, y0_year, k1_year, k2_year = popt_year
    
    ax.plot(x_year[x_year <= x0_year], k1_year*(x_year[x_year <= x0_year]-x0_year)+y0_year, 'r--')
    ax.plot(x_year[x_year >= x0_year], k2_year*(x_year[x_year >= x0_year]-x0_year)+y0_year, 'b--')
    ax.scatter(x0_year, y0_year, c='k')

4. 优化提示 #

• 若数据噪声大，先做滚动平滑处理（如your_df['y_col'].rolling(window=3).mean()）再拟合
• 拟合时指定bounds参数限制斜率范围，比如限制平台段斜率接近0：bounds=([x.min(), y.min(), -np.inf, -0.1], [x.max(), y.max(), 0, 0.1])

内容的提问来源于stack exchange，提问作者GPhil