Newer
Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
# TO DO: THERE IS A BIAS FOR THE LONG LENGHTSCALES
# I GUESS MUST MULTIPLY BY (N-1)/N where N is
# length / scale
#' Haar fluctuation spectrum
#'
#' @description Haar fluctuation spectrum after Lovejoy
#' still needs to be documented
#' currently requires n=2^i data points
#' @param x the input series (typically a numeric)
#' @param q the fluctuation order
#' @examples
#' x = rnorm(2048)
#' xi1 = haar(x, q=1)
#' xi2 = haar(x, q=2)
haar <- function(x, q=2, discarded_scales = 4)
{
n = length(x)
nsteps=floor(log(n)/log(2))
nmax = 2^nsteps
DT = sapply (seq(nsteps), function(i)
{
m = 2^i
nn = nmax / m
M = matrix(x, nmax/m, m, byrow=TRUE)
M1 = M[,1:(m/2), drop=FALSE]
M2 = M[,(m/2+1):m, drop=FALSE]
tM1 = apply(M1, 1, mean)
tM2 = apply(M2, 1, mean)
abs(tM2 - tM1) * ( 2 * nn ) / (2 * nn - 1)
})
x= 2^seq(nsteps)
y= sapply(DT, function(x) mean(x^q))
H = data.frame(x = x, y=y, logx=log(x), logy=log(y))
kept_scales <- seq(max(2,nsteps-discarded_scales))
out = lm( H$logy[kept_scales] ~ H$logx[kept_scales] )$coefficients[2]
attr(out, "scale") = H
attr(out, "q") = q
attr(out,"class") = "fluctuation spectrum"
out
}