从Matalb到Julia的MF-DFA实施

Question

我在Matlab中使用MF-DFA方法，但我需要在Julia中实现它。目标是获得S&第500页。matlab代码如下：

sp500 = readtable('sp500_Nasdaq.csv','PreserveVariableNames', true) ;
spClose = table2array(sp500(:,2))  ;
SP1=cumsum(spClose - mean(spClose)) ;
SP1_ordinary=sqrt(mean(SP1.^2));
X=cumsum(spClose-mean(spClose));
X=transpose(X);
scale=[16,32,64,128,256,512,1024];
q=[-5,-3,-1,0,1,3,5];
m=1;
for ns=1:length(scale),
segments(ns)=floor(length(X)/scale(ns));
for v=1:segments(ns)
Index=( ( ( (v-1)*scale(ns) )+1):(v*scale(ns)));
C = polyfit(Index,X(Index),m) ;
fit=polyval(C,Index);
RMS{ns}(v)=sqrt(mean((X(Index)-fit).^2));
end
for nq=1:length(q),
qRMS{nq,ns}=RMS{ns}.^q(nq);
Fq(nq,ns)=mean(qRMS{nq,ns}).^(1/q(nq));
end
Fq(q==0,ns)=exp(0.5*mean(log(RMS{ns}.^2)));
end

Julia的代码看起来像这样：

using DelimitedFiles, TimeSeries, Plots, DelimitedFiles, Plots, StatsBase
using Polynomials, LinearAlgebra, CSV, DataFrames
sp500 = CSV.read("sp500_Nasdaq.csv", DataFrame) 
sp500_V = values(sp500[:,2]) 
SP1 = cumsum(sp500_V .- mean(sp500_V) ) ; 
SP1_Ord = sqrt(mean(SP1.^2)) ;
X = SP1 ;
X = X';
function polyfit(xVals,yVals)
n = length(xVals)
xBar, yBar = @fastmath mean(xVals), mean(yVals)
sXX, sXY = @fastmath ones(n)'*(xVals.-xBar).^2 , dot(xVals.-xBar,yVals.-yBar)
b1A = @fastmath sXY/sXX
b0A = @fastmath yBar - b1A*xBar
return b0A, b1A
end
scales = [16,32,64,128,256,512,1024];
q = [-5,-3,-1,0,1,3,5] ;
segments = zeros(Int64, (1,length(scales)))
global qRMS = zeros( length(q) ,length(scales)  ) ;
global Fq = zeros( length(q) , length(scales) ) ;
@inbounds for ns = 1:length(scales)
global segments[ns] = Int(floor( length(X)/scales[ns] ) ) ;  
global Index = Array{UnitRange{Int128}}(undef, (segments[ns], length(scales))  ) ;
global ft = zeros(Float64, (segments[ns], length(scales) ) ) ;
global RMS = zeros(Float64, (length(scales) ,segments[ns] ) ) ;
@inbounds  for v=1:segments[ns]
global RMSk = Array{Float64}[] ;
Index =  ( ( (v-1)*scales[ns] ) + 1 ):( v*scales[ns] ) ;
global C = polyfit( Index, X[Index]) ;
global p = Polynomial(C)
ft =p.(Index);
RMS[ns,v] = sqrt(mean((X[Index] .- ft).^2));
push!(RMSk,RMS )
end
@inbounds for nq = 1:length(q)
qRMS[nq,ns] = RMS[ns].^q[nq];
Fq[nq,ns] = mean( qRMS[nq,ns] ).^(1/q[nq] );
end
Fq[findall(x->x==0, q)[1], ns] = exp( 0.5*mean(log.(RMS[ns].^2) ) ) ;
end

问题是，Matlab代码中的阵列RMS是这样的阵列阵列：

RMS =

1×7单元阵列

{1×159 double}    {1×79 double}    {1×39 double}    {1×19 double}    {1×9 double}    {1×4 double}    {1×2 double}

但Julia只返回最后一个数组

RMS
7×2 Matrix{Float64}:
0.0      0.0
0.0      0.0
0.0      0.0
0.0      0.0
0.0      0.0
0.0      0.0
62178.0  18238.2

如何获得与Matlab中相同的输出？如何在Julia中将数组存储到数组中？

Answer 1

这方面的解决方案是使用RMScell=Array｛Float64｝[]等效于Matlab单元阵列

using DelimitedFiles, TimeSeries, Plots, DelimitedFiles, Plots
using Polynomials, LinearAlgebra, CSV, DataFrames, StatsBase
sp500 = CSV.read("sp500_Nasdaq.csv", DataFrame) ;
sp500_V = values(sp500[:,2]) ;
SP1 = cumsum( sp500_V .- mean(sp500_V) ) ; 
SP1_Ord = sqrt( mean(SP1.^2) ) ;
X = SP1 ;
X = X' ;
function polyfit(xVals,yVals)
n = length(xVals)
xBar, yBar = @fastmath mean(xVals), mean(yVals)
sXX, sXY = @fastmath ones(n)'*(xVals.-xBar).^2 , dot(xVals.-  xBar,yVals.-yBar)
b1A = @fastmath sXY/sXX
b0A = @fastmath yBar - b1A*xBar
return b0A, b1A
end
""" Multifractal detrended fluctuation analysis of time series """
scales = [16,32,64,128,256,512,1024];
q = [-5,-3,-1,0,1,3,5] ;
segments = zeros(Int64, (1,length(scales))) ;
global qRMS = zeros( length(q) ,length(scales)  ) ;
global Fq = zeros( length(q) , length(scales) ) ;
global RMScell = Array{Float64}[] ;
global qRMScell =[] ;
global segmentsFq = [] ;
@inbounds for ns = 1:length(scales)
global segments[ns] = Int(floor( length(X)/scales[ns] ) ) ;  
global ft = zeros(Float64, (segments[ns], length(scales) ) ) ;
global RMS = zeros(Float64, segments[ns]);
@inbounds  for v=1:segments[ns]
global Index = ( (v-1)*scales[ns] ) + 1: v*scales[ns] ;
global C = polyfit( Index, X[Index]) ;
global p = Polynomial(C) ;
ft =p.(Index) ;
RMS[v] = sqrt(mean((X[Index] .- ft).^2))  ;
end
l = deepcopy(RMS)
push!(RMScell,l)
global IndexFq = ((ns-1)*length(q) ) + 1 : ns*length(q) ;
push!(segmentsFq, IndexFq) ;
@inbounds for nq = 1:length(q)
l = RMScell[ns].^q[nq]
r = deepcopy(l) ;
push!(qRMScell, r) ;
end
@inbounds for nq = 1: length(scales)
Fq[nq,ns] = mean( qRMScell[segmentsFq[ns]][nq] ).^(1/q[nq] ) ;
end
Fq[findall(x->x==0, q)[1], ns] = exp( 0.5*mean(log.(RMScell[ns].^2) ) ) ;

结束

Hq = zeros( Float64,length(q) ) ;
global qRegLine = Array{Float64}[] ;
for  nq = 1:length(q)
global C = polyfit( log2.(scales),log2.(Fq[nq,:]) ) ;
Hq[nq] = C[2] ;
global p = Polynomial(C) ;
push!( qRegLine, p.( log2.(scales) ) )
end

tq = Hq.*q .- 1 ;
hq = diff(tq)./(q[2]-q[1]) ;
Dq = ( q[1:end-1].*hq ) - tq[1:end-1] ;