平滑python中分段回归载荷-位移曲线的尖角和跳跃

Question

我在平滑模拟软件不太喜欢的一些尖角方面遇到了一个顽固的问题。

我有以下位移/负载/损坏与步骤/时间的关系：

源数据可以在这里找到。

这是导入数据并绘制上图的代码：

df = pd.read_csv("ExampleforStack.txt") # read data
x = df["Displacement"] # get displacement
y = df["Load"] # get load
d = df["Damage"] # get damage
# plot stuff
plt.figure()
plt.subplot(3,1,1)
plt.plot(x)
plt.grid()
plt.ylabel("Displacement")
plt.subplot(3,1,2)
plt.plot(y)
plt.grid()
plt.ylabel("Load")
plt.subplot(3,1,3)
plt.plot(d)
plt.grid()
plt.ylabel("Damage")
plt.xlabel("Step")
plt.gcf().align_ylabels()
plt.tight_layout()

根据位移绘制时，载荷和损坏如下所示：

上图中的断点是：

print(bps)
# [0.005806195310298627, 0.02801208361344569]

我的目标是平滑垂直黑线周围的负载和损坏数据。

到目前为止，我尝试了

statsmodels.api.nonparametric

的lowess，结果看起来非常次优：

上图是 frac 为 0.03 的情况，改变 frac 当然会改变很多，但遗憾的是也不是以理想的方式。

我尝试过的其他东西有高斯回归模型、奇异谱分析、Savitzky-Golay 滤波器、三次样条等...

到目前为止我唯一没有检查的是曲线拟合，我明天可能会检查。

背景资料：

• 位移是DIC分析的结果
• 负载由试验机测量
• 损坏是根据弹性区域中材料的位移、载荷和刚度计算得出的值。

从定性上来说，我希望最终结果如下：

另一个要求是平滑数据的导数也应该平滑且不跳跃。

如果有任何提示可以帮助我解决此任务，我将不胜感激！ ：D

Answer 1

您可以使用（二次）贝塞尔曲线来平滑两个重要点。

在下面的代码中，a、b、c 是定义贝塞尔曲线的三个点的索引。我有点随意地将它们摘下来，以便 b 值与您的断点一致。您可能想要自动化。

import pandas as pd
import matplotlib.pyplot as plt

def Fixit( x, y, d, a, b, c ):
    ''' Quadratic Bezier curve using points with indices a, b, c '''
    N = c - a
    xa, xb, xc = x[a], x[b], x[c]
    ya, yb, yc = y[a], y[b], y[c]
    da, db, dc = d[a], d[b], d[c]
    for i in range( 1, N ):
        t = i / N
        x[a+i] = ( 1 - t ) ** 2 * xa + 2 * t * ( 1 - t ) * xb + t ** 2 * xc
        y[a+i] = ( 1 - t ) ** 2 * ya + 2 * t * ( 1 - t ) * yb + t ** 2 * yc
        d[a+i] = ( 1 - t ) ** 2 * da + 2 * t * ( 1 - t ) * db + t ** 2 * dc




df = pd.read_csv("ExampleforStack.txt") # read data
x = df["Displacement"] # get displacement
y = df["Load"] # get load
d = df["Damage"] # get damage


a, b, c = 210, 240, 250
Fixit( x, y, d, a, b, c )              # Smooth with a Bezier curve

a, b, c = 470, 480, 510
Fixit( x, y, d, a, b, c )              # Smooth with a Bezier curve



# plot stuff
def plotall():
    plt.figure()
    plt.subplot(3,1,1)
    plt.plot(x)
    plt.grid()
    plt.ylabel("Displacement")
    plt.subplot(3,1,2)
    plt.plot(y)
    plt.grid()
    plt.ylabel("Load")
    plt.subplot(3,1,3)
    plt.plot(d)
    plt.grid()
    plt.ylabel("Damage")
    plt.xlabel("Step")
    plt.gcf().align_ylabels()
    plt.tight_layout()
    plt.show()

def plotStress():
    plt.plot(x,y,'b')
    plt.xlabel("Displacement")
    plt.ylabel("Stress")
    plt.show()

def plotDamage():
    plt.plot(x,d,'r')
    plt.xlabel("Displacement")
    plt.ylabel("Damage")
    plt.show()

plotall()
plotStress()
plotDamage()

输出：