размывание работы с Интерполированными функциями в системе Mathematica

Question

Я использую Mathematica 7.

У меня есть интерполированное функцию, вот пример:

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

Я хотел бы вычислить его производный, которая является прямой вперед:

dpressures = D[pressures[x], x] 

Теперь, если вы построить этот Funciton

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

(извините, не знаю, как опубликовать изображение из Mathematica, и у вас нет времени разобраться.) Вы найдете, что это очень шумно. Поэтому я хотел бы сгладить это. Моя первая мысль была использовать свертку, и интегрировать его с гауссовым ядром, что-то вроде следующего:

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.

Answer 1

размывания ищет замкнутую форму для свертки. Вы можете попробовать числовую свертку, начиная с чем-то вроде

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 для получения помех.