cowboy313
一个S语言求解分形高斯噪声的自相似度的程序,请求给予语法上的指点. 本人在网上下载了一个S语言的程序,想转换到matlab代码, 但我没有S语言的基础, 请高手帮忙指点一下.
今天看了一下S语言入门, 觉得程序中第一个函数是不是有问题, 下划线是不是写错了, 应为左箭号,即赋值. 下面贴出前两个函数,如有人帮忙,请留下E-mail, 我将整个程序发到你的邮箱, 谢谢.
###############################################
#
#Splus function and programs for the
#calculation of Whittle's estimator and
#the goodness of fit statistic as defi-
#ned in Beran (1992). The models are
#fractional gaussian noise or fractional
#Arima.
#The data series may be divided into
#subseries for which the parameters are
#fitted separately.
#
###############################################
#
#Functions
#
###############################################
#CetaFGN
#
#Covariance matrix of hat{eta} for fgn
###############################################
CetaFGN_function(eta)
{
M_length(eta)
#size of steps in Riemann sum: 2*pi/m
m_10000
mhalfm_trunc((m-1)/2)
#size of delta for numerical calculation of derivative
delta_0.000000001
#partial derivatives of log f(at each fourier frequency)
lf_matrix(1,ncol=M,nrow=mhalfm)
f0_fspecFGN(eta,m)$fspec
for(j in (1:M))
{
etaj_eta
etaj[j]_etaj[j]+delta
fj_fspecFGN(etaj,m)$fspec
lf[,j]_log(fj/f0)/delta
}
#Calculate D
Djl_matrix(1,ncol=M,nrow=M)
for(j in (1:M))
{for(l in (1:M))
{Djl[j,l]_2*2*pi/m*sum(lf[,j]*lf[,l])
}
}
#Result
drop(matrix(4*pi*solve(Djl),ncol=M,nrow=M,byrow=T))
}
###############################################
#CetaARIMA
#
#Covariance matrix of hat{eta}
#for fractional ARIMA
###############################################
CetaARIMA_function(eta,p,q)
{
M_length(eta)
#size of steps in Riemann sum: 2*pi/m
m_10000
mhalfm_trunc((m-1)/2)
#size of delta for numerical calculation of derivative
delta_0.000000001
#partial derivatives of log f (at each fourier frequency
lf_matrix(1,ncol=M,nrow=mhalfm)
f0_fspecARIMA(eta,p,q,m)$fspec
for(j in (1:M))
{
etaj_eta
etaj[j]_etaj[j]+delta
fj_fspecARIMA(etaj,p,q,m)$fspec
lf[,j]_log(fj/f0)/delta
}
#Calculate D
Djl_matrix(1,ncol=M,nrow=M)
for(j in (1:M))
{for(l in (1:M))
{Djl[j,l]_2*2*pi/m*sum(lf[,j]*lf[,l])
}
}
#Result
drop(matrix(4*pi*solve(Djl),ncol=M,nrow=M,byrow=T))
}