Extrair pontos de dados da média móvel?

15

É possível extrair pontos de dados da movimentação de dados médios?

Em outras palavras, se um conjunto de dados tiver apenas médias móveis simples dos 30 pontos anteriores, é possível extrair os pontos de dados originais?

Se sim, como?


fonte
11
A resposta é sim qualificada, mas o procedimento exato depende de como o segmento inicial de dados é tratado. Se for simplesmente descartado, você perderá efetivamente 15 dados, deixando-o com um sistema indeterminado de equações lineares. O resultado é que existem muitas respostas válidas em geral, mas você ainda pode fazer algum progresso se (a) janelas mais curtas (ou algum procedimento desse tipo) forem usadas para as 15 médias móveis iniciais ou (b) você pode especificar restrições adicionais sobre a solução (cerca de 15 dimensões de restrições ...). Em que situação você está?
whuber
@whuber Muito obrigado por olhar! Eu tenho 2.000 pontos. O primeiro ponto MA é provavelmente uma média dos 30 primeiros pontos originais. A precisão é o segundo resultado geralmente correto, mais especificamente boas suposições nos pontos mais "recentes". Você pode recomendar um método relativamente simples? Desde já, obrigado!
11
(se você demorar mais de cinco minutos para escrever um comentário ...). O que eu queria escrever é que você possa pensar na média como uma multiplicação de matrizes. As linhas no meio terão 1/30 * [1 1 1 ...] antes da diagonal. A questão é: como você lida com pontos nas bordas do seu vetor para tornar a matriz invertível. Você pode fazer isso assumindo que eles resultam da média de menos elementos ou que pensa em outras restrições. Observe que, embora uma inversão de matriz seja uma maneira fácil de entendê-la, ela não é a mais eficiente. Você provavelmente deseja usar uma FFT para fazer isso.
fabee

Respostas:

4

+1 na resposta da fabee, que está completa. Apenas uma nota para traduzi-lo para R, com base nos pacotes que encontrei para executar as operações em mãos. No meu caso, eu tinha dados que são previsões de temperatura da NOAA em uma base de três meses: jan-fev-mar, fev-mar-abr, mar-abr-maio, etc., e queria dividi-los em (aproximado) mensais, assumindo que a temperatura de cada período de três meses seja essencialmente uma média.

library (Matrix)
library (matrixcalc)

# Feb-Mar-Apr through Nov-Dec-Jan temperature forecasts:

qtemps <- c(46.0, 56.4, 65.8, 73.4, 77.4, 76.2, 69.5, 60.1, 49.5, 41.2)

# Thus I need a 10x12 matrix, which is a band matrix but with the first
# and last rows removed so that each row contains 3 1's, for three months.
# Yeah, the as.matrix and all is a bit obfuscated, but the results of
# band are not what svd.inverse wants.

a <- as.matrix (band (matrix (1, nrow=12, ncol=12), -1, 1)[-c(1, 12),])
ai <- svd.inverse (a)

mtemps <- t(qtemps) %*% t(ai) * 3

O que funciona muito bem para mim. Obrigado @fabee.

EDIT: OK, traduzindo meu R para Python, recebo:

from numpy import *
from numpy.linalg import *

qtemps = transpose ([[46.0, 56.4, 65.8, 73.4, 77.4, 76.2, 69.5, 60.1, 49.5, 41.2]])

a = tril (ones ((12, 12)), 2) - tril (ones ((12, 12)), -1)
a = a[0:10,:]

ai = pinv (a)

mtemps = dot (ai, qtemps) * 3

(O que demorou muito mais para depurar do que a versão R. Primeiro, porque eu não estou tão familiarizado com Python quanto com R, mas também porque R é muito mais utilizável interativamente.)

Wayne
fonte
@Gracchus: Desculpe, não é um cara de C ++, mas você pode encontrar o que precisa na biblioteca de álgebra linear do Armadillo C ++ ( arma.sourceforge.net ), que também está disponível no R pelo pacote RcppArmadillo.
Wayne
OK, veja se funciona para você. Se assim for, você poderia escolher minha resposta ;-)
Wayne
As melhores práticas da FYI no Python são fazer importações absolutas: python.org/dev/peps/pep-0008/#imports, o que facilita muito a leitura do código de outras pessoas, porque você realmente sabe de onde vêm as funções em vez de precisar procure cada um que você não conhece. Gostaria que fosse padrão no R fazer o mesmo. Ter que pesquisar cada pequena funções em alguém do código realmente mói meu artes ...
wordsforthewise
Além disso, os notebooks Jupyter para interatividade Python ou IPython.
wordsforthewise
17

Eu tento colocar o que o whuber disse em resposta. Digamos que você tenha um vetor grande com n = 2000 entradas. Se você calcular uma média móvel com uma janela de comprimento = 30 , poderá escrever isso como uma multiplicação da matriz vetorial y = A x do vetor x com a matrizxn=2000=30y=Axx

A=130(1...10...001...10...0...1...100...01...1)

que possui que são deslocados à medida que você avança pelas linhas até que os 30 atinjam o final da matriz. Aqui, o vetor médio y tem dimensões de 1970. A matriz tem 1970 linhas e3030y1970 colunas. Portanto, não é invertível.2000

x1,...,x2000y1y2

x1,...,xnxyx

A3030AA

AAz=AyxyAz (consulte a Wikipedia ).

2000x através do pseudo-inverso.

reconstruction of original signal from moving average using the pseudoinverse

Muitos programas numéricos oferecem pseudo-inversos (por exemplo, Matlab, numpy em python, etc.).

Aqui seria o código python para gerar os sinais do meu exemplo:

from numpy import *
from numpy.linalg import *
from matplotlib.pyplot import *
# get A and its inverse     
A = (tril(ones((2000,2000)),-1) - tril(ones((2000,2000)),-31))/30.
A = A[30:,:]
pA = pinv(A) #pseudo inverse

# get x
x = random.randn(2000) + 5
y = dot(A,x)

# reconstruct
x2 = dot(pA,y)

plot(x,label='original x')
plot(y,label='averaged x')
plot(x2,label='reconstructed x')
legend()
show()

Espero que ajude.

fabee
fonte
Essa é uma ótima resposta, mas acho que você se enganou quando disse que "minimiza a distância quadrática entre y e Az". De fato, y e Az são a mesma coisa. O que é minimizado é a norma de z, que funciona bem para os sinais do mundo real que tentei, mas não é tão boa se o seu sinal original tiver muitos valores discrepantes.
Gdelfino
Não tenho certeza se eu sigo. y e Ax são a mesma coisa, mas não y e Az É verdade que também minimiza a norma de z. Também não vejo por que não funciona nos meus exemplos. A linha azul e a vermelha combinam muito bem. Estou faltando alguma coisa no seu comentário?
Fab28
y é a média móvel calculada a partir do sinal original x multiplicando por A. Este procedimento nos fornece um sinal z que tem a mesma média móvel y. Portanto, y = Az Portanto, apenas a norma de z é minimizada. Se o sinal original tiver um valor padrão alto, o procedimento não dará bons resultados. Um sinal de exemplo com grande valor norma é abaixo:
gdelfino
{42,8, -33,7, 13,2, -45,6, 10,2, 35,8, -41,4, 20.253, 43,3429, -33.2735, 13.6135, -45.1067, 10.6346, 36.1352, -40.9703, 20.6616, 43.6796, -32.8966, 14.0406, -44.7001, 10.9988 , 36.4675, -40.7277, 20.8823, 43.7878, -32.7415, 13.9951, -44.7947, 11.044, 36.3873, -40.7117, 20.7505, 43,8204, -32.9399, 13.9129, -44.9549, 10.8703, 36.1559, -40.8894, 43.484 , 13,5468, -45,2374, 10,3787, 35,8235, -41,5161, 19,9717, 43,0658, -33,7125, 13,0321}
gdelfino
Por favor, use um tamanho de janela de 8 para o sinal acima. Dessa forma, o sinal filtrado tem uma forma muito diferente do sinal original.
Gdelfino