Unclear periodicity output with daily temperature dataset using periodogram
2 次查看(过去 30 天)
显示 更早的评论
Dear all, I am struggling to find a meaningful periodicity of my dataset. Your help would be much appreciated.
I have a dataset (Temp) containing 10 years of temperatures above background. The sampling frequency is one datapoint per day.
Considering that Temp contains the delta T above background, most of the seasonality is smoothened, altought still visible.
As I have several datasets to analyse (where the periodicity might vary), I need to find a way to establish the periodicity of a given dataset.
I tried running the code below (based on https://it.mathworks.com/help/signal/ug/find-periodicity-using-frequency-analysis.html), but I get an unexpected and unclear result.
Time= %MyDailyTempData; 4032x1 double
Temp= %MyDailyTimeData; 4032x1 datetime
fs = 1; % Sampling frequency (1 acquisition per day)
t = (0:length(Temp) - 1)/fs; % Time vector
tempnorm = Temp - mean(Temp); % Subtract the mean to concentrate on temperature fluctuations (based on the link above)
[pxx,f] = periodogram(tempnorm,[],[],fs);
% To find the peak (periodicity of my data)
j=max(pxx);
jj=find(pxx==j);
P=f(jj);
% Plot
figure
subplot(2,1,1)
plot(Time,Temp)
title('Temperature - Mean')
xlabel('Time (Years)')
ylabel('Delta Temperature ( {}^\circC )')
axis tight
subplot(2,1,2)
plot(f,pxx)
grid
title('Periodicity',P)
xlabel('Frequency')
ylabel('Magnitude')
From the above, I obtain:
Now, I would expect my data to have a periodicity of 1 year (by looking at the graph). Additionally, I don't get what 0.1875 represents.
Could you please help me to figure out what I am doing wrong and how to fix it to obtain meaningful results? Ideally, I would like to obtain a bottom graph with the scale in months, or in years. With the pick (in this occasion) heaving a meaningful value of 12 (if months), or 1 (if years).
Thanks in advance for your help!
0 个评论
采纳的回答
Star Strider
2023-2-4
I have no idea what temperature is being measured. I would expect environmental temperatures to have a significantly wider range than about ±4°C.
If the sampling frequency is 1/day, the sampling interval is 1 day, and that means the frequency axis is in cycles/day. (The frequency axis extends from 0 cycles/day to the Nyquist frequency of 0.5 cycles/day.)
The frequency of that peak is 0.1875 cycles/day. The inverse of that is a period of about 5.3 days/cycle. I have no idea how to interpret that, since I know nothing about the location of the thermometer, or what else might be occurring there.
The small peak at the left of the Fourier transform plot (that appears to be about 0.005 cycles/day) might make more sense. It corresponds to a period of about 200 days. (A period of 365 days corresponds to a frequency of about 0.0027 cycles/day, so that might actually be the small peak. I cannot determine exactly what its frequency is.)
2 个评论
Star Strider
2023-2-4
The inversion is straightforward.
Specifically:
so:
They are simply the inverses of each other.
And the inverses of 0.00024 and 0.00268 are respectively:
format shortG
periods = 1./[0.00024 0.00268]
in_years = periods/365.25
So the lowest frequency is about 11 years, approximately corresponding to the sunspot cycle, and the next highest frequency is just about 1 year. (I still do not understand the peak frequency that corresponds to a period of about 5.3 days. Perhaps it has something to do with the satellite or its orbit.)
I do not usually use the periodogram function so I have limited experience with it. From the documentation: ‘The units of the PSD estimate are in squared magnitude units of the time series data per unit frequency.’ This (to me) makes the results a bit more difficult to interpret than the results of the fft or pspectrum functions, that I generally prefer, where the units are simply in terms of the magnitude (or, if squared, power) of the signal at specific frequencies, rather than in power per frequency unit.
更多回答(1 个)
Image Analyst
2023-2-4
I agree with everything Star said - I was going to say the same things. But what is "background"? Does that change seasonally? If so, maybe "delta T" is fairly constant seasonally????
I don't really have suggestions for what a 5 day periodicity in the data might mean. Now 7 days I could possibly see because that coincides with typical work weeks and the temperature may vary periodically weekly as people stay home for the weekend.
另请参阅
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!