如何使用 Matlab 找到噪声数据系列的拐点？

Question

给定以下曲线

我想确定曲线开始真正增加的 x 数据点的索引（在本例中，大约是 x=15）。

虽然我知道可以使用导数来确定拐点，但请注意，数据是嘈杂的，我不确定这种方法能否让我清楚地识别“真正的拐点”（在这种情况下 x=15）。

我想知道一个更简单的方法是否可行，例如

• 找到 4 个数据点，其中 x1
• x1 的返回索引

您对如何完成此任务有任何建议吗？

以上曲线的样本数据

 index       SQMean   
_____    ____________

'0'      '139.428574'
'1'      '133.298706'
'2'      '135.961044'
'3'      '143.688309'
'4'      '133.298706'
'5'      '133.181824'
'6'      '134.896103'
'7'      '146.415588'
'8'      '142.324677'
'9'      '128.168839'
'10'     '146.116882'
'11'     '146.766235'
'12'     '134.675323'
'13'     '138.610382'
'14'     '140.558441'
'15'     '128.662338'
'16'     '138.480515'
'17'     '153.610382'
'18'     '156.207794'
'19'     '183.428574'
'20'     '220.324677'
'21'     '224.324677'
'22'     '230.415588'
'23'     '226.766235'
'24'     '223.935059'
'25'     '229.922073'
'26'     '234.389618'
'27'     '235.493500'
'28'     '225.727280'
'29'     '241.623383'
'30'     '225.805191'
'31'     '240.896103'
'32'     '224.090912'
'33'     '230.467529'
'34'     '248.285721'
'35'     '233.779221'
'36'     '225.532471'
'37'     '247.337662'
'38'     '233.000000'
'39'     '241.740265'
'40'     '235.688309'
'41'     '238.662338'
'42'     '236.636368'
'43'     '236.025970'
'44'     '234.818176'
'45'     '240.974030'
'46'     '251.350647'
'47'     '241.857147'
'48'     '242.623383'
'49'     '245.714279'
'50'     '250.701294'
'51'     '229.415588'
'52'     '236.909088'
'53'     '243.779221'
'54'     '244.532471'
'55'     '241.493500'
'56'     '245.480515'
'57'     '244.324677'
'58'     '244.025970'
'59'     '231.987015'
'60'     '238.740265'
'61'     '239.532471'
'62'     '232.363632'
'63'     '242.454544'
'64'     '243.831161'
'65'     '229.688309'
'66'     '239.493500'
'67'     '247.324677'
'68'     '245.324677'
'69'     '244.662338'
'70'     '238.610382'
'71'     '243.324677'
'72'     '234.584412'
'73'     '235.181824'
'74'     '228.974030'
'75'     '228.246750'
'76'     '230.519485'
'77'     '231.441559'
'78'     '236.324677'
'79'     '229.935059'
'80'     '238.701294'
'81'     '236.441559'
'82'     '244.350647'
'83'     '233.714279'
'84'     '243.753250'

Answer 1

如果这是一次性估计，您可以做的一件事是使用曲线拟合工具箱中的曲线拟合工具。这是一个示例，我将分段线性函数拟合到您的数据：

（点击图片查看大图）

函数形式

a * (x < b) + c * (x > d) + ((x - b) / (d - b) * (c - a) + a) * (x >= b) * (x <= d)

它说：x < b 有一个常数部分，其值为 a，x > d 的另一个常数部分，其值为 c，以及连接它们的线性斜坡。

很难适应这样的功能，并且只有在您提供合适的起始估计时才能正常工作（请参阅屏幕截图中的小窗口）。因此，这不是一种使流程自动化的方法，而只是为了获得改进的估计。

在这种情况下，从 b = 15 的初始估计值开始，拟合提供了改进的估计值 b = 16.58，其 95%-CI 为 [15.96, 17.2]，这表明索引 0 到 16 属于初始常数部分。

曲线拟合工具还可以根据您的 GUI 规范生成代码。在这种情况下，结果是：

[xData, yData] = prepareCurveData( index, SQMean );

% Set up fittype and options.
ft = fittype( 'a * (x < b) + c * (x > d) + ((x - b) / (d - b) * (c - a) + a) * (x >= b) * (x <= d)', 'independent', 'x', 'dependent', 'y' );
opts = fitoptions( ft );
opts.Display = 'Off';
opts.Lower = [-Inf -Inf -Inf -Inf];
opts.StartPoint = [140 15 230 20];
opts.Upper = [Inf Inf Inf Inf];

% Fit model to data.
[fitresult, gof] = fit( xData, yData, ft, opts );

% Plot fit with data.
figure( 'Name', 'untitled fit 1' );
h = plot( fitresult, xData, yData );
legend( h, 'SQMean vs. index', 'untitled fit 1', 'Location', 'NorthEast' );
% Label axes
xlabel( 'index' );
ylabel( 'SQMean' );
grid on