Seasonal decomposition of a daily time series
21 views (last 30 days)
Show older comments
I have a time series that contains daily variation from 2015 to 2023. I want to see the trend, seasonal components and unmodelled part. Kindly help. The time series has veen attached
0 Comments
Answers (2)
Steven Lord
on 2 Aug 2024 at 14:57
Take a look at the trenddecomp function (introduced in release R2021b) and/or the Find and Remove Trends Live Editor Task (introduced in release R2019b.)
3 Comments
Steven Lord
on 2 Aug 2024 at 16:24
Which release are you using?
You could try using the detrend function as part of writing your own seasonal decomposition code. That function was introduced sometime prior to release R2006a.
Star Strider
on 2 Aug 2024 at 17:12
Your data are not regularly-sampled, so it is necessary to use the funciton on them first. After that, calculate the Fourier transform of the data to see the relevant peaks. From there, you can use the lowpass or bandpass functions as appropriate to get the relevant information from your data.
Try this —
format longG
A1 = readmatrix('time_series.txt.txt')
Ts = mean(diff(A1(:,1)))
Fs = 1/Ts
Tsd = std(diff(A1(:,1)))
[A2(:,2), A2(:,1)] = resample(A1(:,2), A1(:,1), Fs);
A2
figure
plot(A1(:,1), A1(:,2))
grid
xlabel('Time (years)')
ylabel('Amplitude')
title('Original')
figure
plot(A2(:,1), A2(:,2))
grid
xlabel('Time (years)')
ylabel('Amplitude')
title('Resampled')
[FTs1,Fv] = FFT1(A2(:,2), A2(:,1));
[pks,locsidx] = findpeaks(abs(FTs1)*2, 'MinPeakProminence',2E-3)
figure
plot(Fv, abs(FTs1)*2)
hold on
plot(Fv(locsidx), pks, 'vr')
hold off
grid
xlim([0 10])
xlabel('Frequency (Years^{-1})')
ylabel('Magnitude (Units)')
text(Fv(locsidx), pks, compose('\\leftarrow Frequency = %.3f years^{-1} = (1/%.3f years)', [Fv(locsidx); 1./Fv(locsidx)].') )
function [FTs1,Fv] = FFT1(s,t)
% Arguments:
% s: Signal Vector Or Matrix
% t: Associated Time Vector
t = t(:);
L = numel(t);
if size(s,2) == L
s = s.';
end
Fs = 1/mean(diff(t));
Fn = Fs/2;
NFFT = 2^nextpow2(L);
FTs = fft((s - mean(s)) .* hann(L).*ones(1,size(s,2)), NFFT)/sum(hann(L));
Fv = Fs*(0:(NFFT/2))/NFFT;
% Fv = linspace(0, 1, NFFT/2+1)*Fn;
Iv = 1:numel(Fv);
FTs1 = FTs(Iv,:);
end
Experiment to get the result you want. I will help with the filtering if necessary.
.
0 Comments
See Also
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!