如何在 R 中拟合平滑滞后？答案

【问题标题】：How can I fit a smooth hysteresis in R?如何在 R 中拟合平滑滞后？
【发布时间】：2018-09-11 17:26:03
【问题描述】：

我有一个适合滞后的测量值。出于可视化目的，我想绘制一条近似滞后的线来帮助解释这种模式。

我使用下面的代码在下图中创建了一个示例。

我想要一个类似于绿色曲线的输出 - 但是我没有直接可用的数据，我不在乎它是否是尖的。

但是大多数平滑函数，例如我用蓝色绘制的smooth.spline - 不允许循环。我能找到的最接近的来自bezier 库 - 以红色绘制。在这里不太明显，但它会产生一个循环，但它很不适合（并且会给出一些警告并且需要相当长的时间）。

你能推荐一个方法吗？

set.seed(12345)
up <- seq(0,1,length.out=100)^3
down <- sqrt(seq(1,0,length.out=100))
x <- c(seq(0,1,length.out=length(up)),
       seq(1,0, length.out=length(down)))

data <- data.frame(x=x, y=c(up,down),
                   measuredx=x + rnorm(length(x))*0.01,
                   measuredy=c(up,down) + rnorm(length(up)+length(down))*0.03)


with(data,plot(measuredx,measuredy, type = "p"))
with(data,lines(x,y, col='green'))

sp <- with(data,smooth.spline(measuredx, measuredy))
with(sp, lines(x,y, col="blue"))


library(bezier)
bf <- bezierCurveFit(as.matrix(data[,c(1,3)]))
lines(bezier(t=seq(0, 1, length=500), p=bf$p), col="red", cex=0.25)

更新

事实证明我的实际问题略有不同，我提出另一个问题以反映我在问题中的实际问题：How to fit a smooth hysteresis in a poorly distributed data set?

【问题讨论】：

标签： r spline smoothing

【解决方案1】：

set.seed(12345)
up <- seq(0,1,length.out=100)^3
down <- sqrt(seq(1,0,length.out=100))
x <- c(seq(0,1,length.out=length(up)),
       seq(1,0, length.out=length(down)))

data <- data.frame(x=x, y=c(up,down),
                   measuredx=x + rnorm(length(x))*0.01,
                   measuredy=c(up,down) + rnorm(length(up)+length(down))*0.03)

不要直接在data$measuredx 上平滑data$measuredy，而是进行两次单独的平滑处理，分别针对时间戳变量进行平滑处理。然后结合两个平滑的拟合值。这是平滑闭合曲线或循环的一般方法。（另见问答：Smoothing Continuous 2D Points）

t <- seq_len(nrow(data) + 1)
xs <- smooth.spline(t, c(data$measuredx, data$measuredx[1]))$y
ys <- smooth.spline(t, c(data$measuredy, data$measuredy[1]))$y
with(data, plot(measuredx, measuredy))
lines(xs, ys)

c(data$measuredx, data$measuredx[1])例如只是为了保证向量中的最后一个值与第一个一致，从而完成一个循环。

曲线在左下角并不是真正闭合的，因为smooth.spline是在做平滑而不是插值，所以即使我们保证数据向量完成一个循环，拟合的也不一定是闭合的。一个实用的解决方法是使用加权回归，在这个位置施加很大的权重以使其闭合。

t <- seq_len(nrow(data) + 1)

w <- rep(1, length(t))  ## initially identical weight everywhere
w[c(1, length(w))] <- 100000  ## give heavy weight

xs <- smooth.spline(t, c(data$measuredx, data$measuredx[1]), w)$y
ys <- smooth.spline(t, c(data$measuredy, data$measuredy[1]), w)$y
with(data, plot(measuredx, measuredy), col = 8)
lines(xs, ys, lwd = 2)

【讨论】：