function [F alpha pred] = detrended_fluctuation_analysis(X,nsamps,overlap)
% [F alpha pred] = detrended_fluctuation_analysis(X,nsamps,overlap)
%
% Inputs:
% X - vector of power over time
% nsamps - vector of sample numbers defining the time segments used for the
% analysis; segment lengthes should be logarithmic
% overlap - expressed as proportion, e.g., 0 - no overlap, 0.5 - 50% overlap
%
% Outputs:
% F - dfa fluctuation parameter
% alpha - slope of linear fit to log10(dfa) as a function log10(nsamps);
% alpha ~0.5 --> uncorrelated signal (white noise)
% alpha >0.5 --> temporal correlations over the range of tau
% pred - predicted dfa values from linear fit
%
% Steps:
% (1) mean subtraction, (2) cummulative sum, (3) detrend segments, (4) rms
% of detrended sgement, (5) fitting a linear function to dfa
%
% References:
% Linkenkaer-Hansen K et al. (2001) Long-Range Temporal Correlations and
% Scaling Behavior in Human Brain Oscillations. J Neurosci 21:1370–1377.
% Smit DJA et al. (2013) Long-Range Temporal Correlations in Resting-State Alpha
% Oscillations Predict Human Timing-Error Dynamics. J Neurosci 33:11212–11220.
%
% Description:
% The script calculates long-range temporal correlations of oscillatory
% activity. It relates to the slope of the 1/f scale-free properties of
% power spectra. It is an autocorrelation analysis for non-stationary data.
% ---------
%
% Copyright (C) 2015, B. Herrmann
% This program is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% This program is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with this program. If not, see .
%
% ----------------------------------------------------------------------
% B. Herrmann, Email: herrmann.b@gmail.com, 2015-10-23
% subtract mean from amplitude envelope
X = X(:);
X = X - mean(X);
nsamps = nsamps(:);
% calculate cummulative sum
CX = cumsum(X);
% get shifting vector
shift = round(nsamps*(1-overlap));
F = zeros([length(nsamps) 1]);
for ii = 1 : length(nsamps)
nn = 1:nsamps(ii);
counter = 0;
while nn(end) < length(CX)
% detrend segment
tmpX = detrend(CX(nn));
% calculate RMS
F(ii) = F(ii) + sqrt(sum(tmpX.^2)/length(tmpX));
nn = nn + shift(ii);
counter = counter + 1;
end
F(ii) = F(ii)/counter;
end
coef = polyfit(log10(nsamps),log10(F),1); % linear fit
pred = polyval(coef,log10(nsamps)); % predicted dfa values
alpha = coef(1); % slope