parestimate: Estimate periods from (set of) eigenvectors

Description

Function to estimate the parameters (frequencies and rates) given a set of SSA eigenvectors.

Usage

# S3 method for 1d.ssa
parestimate(x, groups, method = c("esprit", "pairs"),
            subspace = c("column", "row"),
            normalize.roots = NULL,
            dimensions = NULL,
            solve.method = c("ls", "tls"),
            …,
            drop = TRUE)
# S3 method for toeplitz.ssa
parestimate(x, groups, method = c("esprit", "pairs"),
            subspace = c("column", "row"),
            normalize.roots = NULL,
            dimensions = NULL,
            solve.method = c("ls", "tls"),
            …,
            drop = TRUE)
# S3 method for mssa
parestimate(x, groups, method = c("esprit", "pairs"),
            subspace = c("column", "row"),
            normalize.roots = NULL,
            dimensions = NULL,
            solve.method = c("ls", "tls"),
            …,
            drop = TRUE)
# S3 method for cssa
parestimate(x, groups, method = c("esprit", "pairs"),
            subspace = c("column", "row"),
            normalize.roots = NULL,
            dimensions = NULL,
            solve.method = c("ls", "tls"),
            …,
            drop = TRUE)
# S3 method for nd.ssa
parestimate(x, groups,
            method = c("esprit"),
            subspace = c("column", "row"),
            normalize.roots = NULL,
            dimensions = NULL,
            solve.method = c("ls", "tls"),
            pairing.method = c("diag", "memp"),
            beta = 8,
            …,
            drop = TRUE)

Arguments

SSA object

groups

list of indices of eigenvectors to estimate from

…

further arguments passed to 'decompose' routine, if necessary

drop

logical, if 'TRUE' then the result is coerced to lowest dimension, when possible (length of groups is one)

dimensions

a vector of dimension indices to perform ESPRIT along. 'NULL' means all dimensions.

method

For 1D-SSA, Toeplitz SSA, and MSSA: parameter estimation method, 'esprit' for 1D-ESPRIT, 'pairs' for rough estimation based on pair of eigenvectors. For nD-SSA: parameter estimation method. For now only 'esprit' is supported.

solve.method

approximate matrix equation solving method, 'ls' for least-squares, 'tls' for total-least-squares.

pairing.method

method for esprit roots pairing, 'diag' for `2D-ESPRIT diagonalization', 'memp' for ``MEMP with an improved pairing step'

subspace

which subspace will be used for parameter estimation

normalize.roots

logical vector or 'NULL', force signal roots to lie on unit circle. 'NULL' means automatic selection: normalize iff circular topology OR Toeplitz SSA used

beta

In nD-ESPRIT, coefficient(s) in convex linear combination of shifted matrices. The length of beta should be ndim - 1, where ndim is the number of independent dimensions. If only one value is passed, it is expanded to a geometric progression.

Value

For 1D-SSA (and Toeplitz), a list of objects of S3-class `fdimpars.1d'. Each object is a list with 5 components:

roots: complex roots of minimal LRR characteristic polynomial
periods: periods of dumped sinusoids
frequencies: frequencies of dumped sinusoids
moduli: moduli of roots
rates: rates of exponential trend (rates == log(moduli))

For 'method' = 'pairs' all moduli are set equal to 1 and all rates equal to 0.

For nD-SSA, a list of objects of S3-class `fdimpars.nd'. Each object is named list of n `fdimpars.1d' objects, each for corresponding spatial coordinate.

In all cases elements of the list have the same names as elements of groups. If group is unnamed, corresponding component gets name `Fn', where `n' is its index in groups list.

If 'drop = TRUE' and length of 'groups' is one, then corresponding list of estimated parameters is returned.

Details

The time series is assumed to satisfy the model $x_{n} = \sum_{k} C_{k} μ_{k}^{n}$ for complex $μ_{k}$ or, alternatively, $x_{n} = \sum_{k} A_{k} ρ_{k}^{n} \sin (2 π ω_{k} n + ϕ_{k}) .$ The return value are the estimated moduli and arguments of complex $μ_{k}$ , more precisely, $ρ_{k}$ ('moduli') and $T_{k} = 1 / ω_{k}$ ('periods').

For images, the model $x_{i j} = \sum_{k} C_{k} λ_{k}^{i} μ_{k}^{j}$ is considered.

Also `print' and `plot' methods are implemented for classes `fdimpars.1d' and `fdimpars.nd'.

References

Golyandina, N., Zhigljavsky, A. (2013): Singular Spectrum Analysis for time series. Springer Briefs in Statistics. Springer.

Roy, R., Kailath, T., (1989): ESPRIT: estimation of signal parameters via rotational invariance techniques. IEEE Trans. Acoust. 37, 984--995.

Rouquette, S., Najim, M. (2001): Estimation of frequencies and damping factors by two- dimensional esprit type methods. IEEE Transactions on Signal Processing 49(1), 237--245.

Wang, Y., Chan, J-W., Liu, Zh. (2005): Comments on ``estimation of frequencies and damping factors by two-dimensional esprit type methods''. IEEE Transactions on Signal Processing 53(8), 3348--3349.

Examples

Run this code

# NOT RUN {
# Decompose 'co2' series with default parameters
s <- ssa(co2, neig = 20)
# Estimate the periods from 2nd and 3rd eigenvectors using 'pairs' method
print(parestimate(s, groups = list(c(2, 3)), method = "pairs"))
# Estimate the peroids from 2nd, 3rd, 5th and 6th eigenvectors using ESPRIT
pe <- parestimate(s, groups = list(c(2, 3, 5, 6)), method = "esprit")
print(pe)
plot(pe)

# }
# NOT RUN {
# Artificial image for 2D SSA
mx <- outer(1:50, 1:50,
            function(i, j) sin(2*pi * i/17) * cos(2*pi * j/7) + exp(i/25 - j/20)) +
      rnorm(50^2, sd = 0.1)
# Decompose 'mx' with default parameters
s <- ssa(mx, kind = "2d-ssa")
# Estimate parameters
pe <- parestimate(s, groups = list(1:5))
print(pe)
plot(pe, col = c("green", "red", "blue"))

# Real example: Mars photo
data(Mars)
# Decompose only Mars image (without background)
s <- ssa(Mars, mask = Mars != 0, wmask = circle(50), kind = "2d-ssa")
# Reconstruct and plot texture pattern
plot(reconstruct(s, groups = list(c(13,14, 17, 18))))
# Estimate pattern parameters
pe <- parestimate(s, groups = list(c(13,14, 17, 18)))
print(pe)
plot(pe, col = c("green", "red", "blue", "black"))
# }

Run the code above in your browser using DataLab

Last chance! 50% off unlimited learning