在Mathematica中使用卷积处理插值函数

Question

我正在使用Mathematica 7。

我有一个插值函数，以下是一个示例：

pressures = 
  WeatherData["Chicago", "Pressure", {2010, 8}] // 
     DeleteCases[#, {_, _Missing}] & // 
    Map[{AbsoluteTime[#[[1]]], #[[2]]} &, #] & // Interpolation;

我希望能够计算它的导数，这很简单：

dpressures = D[pressures[x], x]

现在，如果您绘制此函数

Plot[3600*dpressures, {x, AbsoluteTime[{2010, 8, 2}], AbsoluteTime[{2010, 8, 30}]}]

抱歉，我不知道如何在Mathematica中从内部发布图像，并且没有时间弄清楚。你会发现它非常嘈杂。因此，我想要平滑它。我的第一个想法是使用Convolve，并将其与高斯核相结合，类似于以下内容：

a = Convolve[PDF[NormalDistribution[0, 5], x], 3600*dpressures, x, y]

返回

360 Sqrt[2/\[Pi]] Convolve[E^(-(x^2/50)), InterpolatingFunction[{{3.48961266 10^9, 3.49228746 10^9}},<>], ][x], x, y]

在我的看法中，这看起来是合理的。不幸的是，我认为我在某个地方犯了一个错误，因为我得到的结果似乎无法被评估。也就是说：

a /. y -> AbsoluteTime[{2010, 8, 2}]

返回结果

360 Sqrt[2/\[Pi]] Convolve[E^(-(x^2/50)), InterpolatingFunction[{{3.48961266 10^9, 3.49228746 10^9}},<>][x], x, 3489696000]]

我希望得到的是一个介于-1和1之间的数字，但这似乎并不是我想要的。

- John

Answer 1

Convolve寻求卷积的闭合形式。您可以尝试数值卷积，从类似以下内容开始：

NConvolve[f_, g_, x_, y_?NumericQ] := 
 NIntegrate[f (g /. x -> y - x), {x, -Infinity, Infinity}]

然而，对于这种嘈杂的非平滑函数，数值积分会遇到困难。(默认设置无法工作，即使仔细选择设置，速度也会很慢。)

我建议您直接操作基础数据而不是插值嘈杂的数据。

您的时间范围限制为:

In[89]:= {lower = Min[First[pressures]], upper = Max[First[pressures]]}    
Out[89]= {3.48961*10^9, 3.49229*10^9}

使用插值法获得每小时样本*：

data = Table[pressures[x], {x, lower, upper, 3600}];

现在进行比较。

ListLinePlot[Differences[data]]

使用5小时时间窗口平滑版本：

ListLinePlot[GaussianFilter[Differences[data], 5]]

• 对于噪声数据，您可能希望使用InterpolationOrder -> 1。

- Andrew Moylan