#' Multiscale Multifractal Analysis
#'
#' Applies the Multiscale Multifractal Analysis (MMA) on time series.
#' @usage MMA(tsx, scale, qminmax, ovlap=0, m=2)
#' @param tsx Univariate time series (must be a vector or a ts object).
#' @param scale Vector of scales.
#' @param qminmax Vector of two values min and max of q-order of the moment.
#' @param ovlap Overlapping parameter (By default ovlap=0: no overlapping).
#' @param m Polynomial order for the detrending (by defaults m=2).
#' @return A matrix with three columns (q-order, scale (s), and the scale exponent).
#'
#' @note The original code of this function is in Matlab, you can find it on the
#' following website \href{https://physionet.org/physiotools/mma/}{Physionet}. See
#' references below.
#'
#'
#'
#' @examples
#'
#' \dontrun{
#' library(MFDFA)
#' library(plotly)
#' library(plot3D)
#'
#' a<-0.6
#' N<-800
#' tsx<-MFsim(N,a)
#' scale=10:100
#' res<-MMA(tsx, scale, qminmax=c(-10,10), ovlap=0, m=2)
#'
#' ## Visualisation 1:
#' S_exponent <- matrix(res[,3], nrow=length(unique(res[,1])), ncol=length(min(scale):(max(scale)/5)))
#' m_scale <- unique(res[,2])
#' q <- unique(res[,1])
#' plot_ly() %>% add_surface(x = ~m_scale, y = ~q,
#' z = ~S_exponent)
#'
#' ## Visualisation 2:
#' image2D(S_exponent, xlab="q", ylab="scale", axes=F)
#' axis(1, seq(0,1,0.1), round(quantile(q, seq(0, 1, 0.1)), 2))
#' axis(2, seq(0,1,0.1), round(quantile(m_scale, seq(0, 1, 0.1)), 2))
#' }
#'
#' @references
#' J. Feder, Fractals, Plenum Press, New York, NY, USA, 1988.
#'
#' J. Gieraltowski, J. J. Zebrowski, and R. Baranowski,
#' Multiscale multifractal analysis of heart rate variability recordings
#' http://dx.doi.org/10.1103/PhysRevE.85.021915
#'
#' Goldberger AL, Amaral LAN, Glass L, Hausdorff JM, Ivanov PCh, Mark RG,
#' Mietus JE, Moody GB, Peng C-K, Stanley HE. PhysioBank, PhysioToolkit,
#' and PhysioNet: Components of a New Research Resource for Complex
#' Physiologic Signals. Circulation 101(23):e215-e220.
#'
#' J. W. Kantelhardt, S. A. Zschiegner, E. Koscielny-Bunde, S. Havlin,
#' A. Bunde, H. Stanley, Multifractal detrended fluctuation analysis of
#' nonstationary time series, Physica A: Statistical Mechanics and its
#' Applications, 316 (1) (2002) 87 – 114.
#'
#' J. Gierałtowski, J. J. Żebrowski, and R. Baranowski, "Multiscale
#' multifractal analysis of heart rate variability recordings with a
#' large number of occurrences of arrhythmia," Phys. Rev. E 85, 021915 (2012)
#'
#'
#' @importFrom numbers mod
#' @export
#'
#'
MMA <- function(tsx, scale, qminmax, ovlap=0, m=2){
qrange <- seq(qminmax[1], qminmax[2], 0.1)
qrange[which(qrange==0)]<-0.0001
prof = cumsum(tsx)
slength = length(prof)
fqs<-c()
for (i in scale){
if (ovlap==1){
vec <- 0:(i-1)
ind <- 1:(slength-i+1)
A<-matrix(rep(vec, each = length(ind)), length(ind), length(vec))
coorxy <- apply(A, 2, FUN=function(x) (x+ind))
} else {
nd<-(slength-mod(slength,i))
nc<-(slength-mod(slength,i))/i
coorxy<-matrix(1:nd, nc, i, byrow = TRUE)
}
segments <- apply(coorxy, 2, FUN=function(x)(prof[x]))
xbs <- 1:i
f2nis <- c()
for (j in 1:nrow(segments)){
seg <- segments[j,]
ft <- poly_fit.val(xbs, seg, m)
fit<-ft$polyval
f2nis <-c(f2nis, mean((seg- fit)^2))
}
for (qq in qrange){
fqs<- c(fqs, qq, i, mean(f2nis^(qq/2))^(1/qq))
}
}
fqs1 <- matrix(fqs, ncol=3, byrow = T)
fqs1 <- matrix(append(fqs1, log(fqs1[,2])), ncol=4)
fqs1 <- matrix(append(fqs1, log(fqs1[,3])), ncol=5)
hqs <- c()
for (i in min(scale):(max(scale)/5)){
for (qq in qrange){
tempfqs <- fqs1[which(fqs1[,1]==qq & fqs1[,2] >= i & fqs1[,2] <= 5*i),]
hft <- poly_fit.val(tempfqs[,4], tempfqs[,5], 1)
hqs <- c(hqs, qq, i, as.numeric(hft$polyfit[1]))
}
}
hqs1 <- matrix(hqs, ncol=3, byrow = T)
colnames(hqs1)<- c("q", "s", "S_Exp")
return(hqs1)
}
poly_fit.val<-function(x,y,n){
formule<-lm(as.formula(paste('y~',paste('I(x^',1:n,')', sep='',collapse='+'))))
res1<-coef(formule)
poly.res<-res1[length(res1):1]
suppressWarnings(res2<-predict(formule, interval = "prediction"))
poly.eva<-res2[,1]
allres<-list(polyfit=round(poly.res,4), polyval=poly.eva)
return(allres)}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.