Unclear periodicity output with daily temperature dataset using periodogram

3 views (last 30 days)
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!

Accepted Answer

Star Strider
Star Strider on 4 Feb 2023
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 Comments
Simone A.
Simone A. on 4 Feb 2023
Dear @Star Strider, Thanks a lot for your reply. This temperature dataset is obtaines from satellite scenes. Each datapoint is the temperature of a given central pixel minus the neighbouring 8 (within a 9x9 matrix).
Spikes in the dataset are likely due to contamination of cold cloudy pixels. The small peaks at the left of the Fourier transform plot are 0.00024 and 0.00268, respectively. See below:
I am looking for changes in trends, thus I need to establish the seasonality for decomposition.
Could you kindly explain, or provide more information, on how to obtain this inversion "The frequency of that peak is 0.1875 cycles/day. The inverse of that is a period of about 5.3 days/cycle.".
Additionally, is there a way to set the graph in days, months or years? Being not an expert in this field, that would make my life easier, and will give me a better chance to understand the dataset.
Kindest regards
Star Strider
Star Strider on 4 Feb 2023
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]
periods = 1×2
1.0e+00 * 4166.7 373.13
in_years = periods/365.25
in_years = 1×2
11.408 1.0216
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.

Sign in to comment.

More Answers (1)

Image Analyst
Image Analyst on 4 Feb 2023
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.
  1 Comment
Simone A.
Simone A. on 4 Feb 2023
Thanks for your reply. The rationale behind the investigation is:
I do expect this given pixel to be hotter that its surrounding. I take the temperature of the surrounding 8 pixels (9x9 matrix excluding the central, "potentially" hotter, pixel). This is my background temperature. I subtract the background averaged temperature from the central (potentially hotter) pixel. This gives me the delta T of my dataset.
Given the approach above, I did not expect much seasonality left, as this should have been mostly removed using this approach. However, this is not a big deal, as long as this can be identified and removed.
The issue is that I have no idea why this period peaks at 5 days, considering that my region of interest is in the middle of nowere, and there are no anthropogenic factors that can affect this dataset. In theory, the only residual fluctuations or periods should be seasonal and harmonical.

Sign in to comment.

Categories

Find more on Time Series in Help Center and File Exchange

Products


Release

R2022b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!