Ajustar dados lineares por partes

18

Qual é uma maneira robusta de ajustar dados lineares, porém barulhentos, por partes?

Estou medindo um sinal, que consiste em vários segmentos quase lineares. Eu gostaria de ajustar automaticamente várias linhas aos dados para detectar as transições.

O conjunto de dados consiste em alguns milhares de pontos, com 1 a 10 segmentos e eu sei o número de segmentos.

Este é um exemplo do que eu gostaria de fazer automaticamente.

insira a descrição da imagem aqui

P3trus
fonte
Eu não acho que essa pergunta possa ser respondida razoavelmente, a menos que você nos diga com que precisão você deseja conhecer os locais dos pontos de interrupção, qual é o seu estimativa de estimativa para o menor comprimento de um segmento linear e quantas amostras existem em um típico região de transição. Se os rótulos dos eixos horizontais na sua figura são números de amostra, então, com duas transições no intervalo de a , a tarefa é mais difícil do que se os segmentos retos tivessem maior duração (em amostras). x [ 0 ]x[5]x[0]
usar o seguinte código
@DilipSarwate Eu atualizei a questão com os requisitos (aliás os xaxis é o campo magnético em tesla)
P3trus
Você pode tentar esta caixa de ferramentas se estiver trabalhando com a caixa de ferramentas para ajuste de curvas
Rhei

Respostas:

12

Eu tentei duas abordagens, ingenuamente (usando apenas 3 segmentos). Certamente haveria métodos mais sofisticados por aí.

    RANSAC, deveria ser um mecanismo de encaixe robusto. É fácil interromper o algoritmo após vários segmentos. No entanto, pode ser difícil impor a continuidade entre segmentos - como parece necessário no seu aplicativo - pelo menos com uma implementação simples. Como prova de conceito, criei uma imagem a partir dos pontos de dados para poder usar o mecanismo RANSAC disponível no , a função de detecção de linha do Mathematica.ImageLines

insira a descrição da imagem aqui

    Ajuste um modelo linear por partes usando um minimizador de uso geral. É fácil impor a continuidade dos segmentos. Curiosamente, o teste de resíduos e outras propriedades pode fornecer informações suficientes para determinar automaticamente o número de segmentos - eu ainda não tentei. É assim que parece no Mathematica:

insira a descrição da imagem aqui

Matthias Odisio
fonte
Parece uma ótima resposta. Obrigado por contribuir.
Jason R
7

Não afirmo que o método a seguir seja robusto, mas pode funcionar para você. Com milhares de pontos e talvez dez ou mais segmentos retos, faça o seguinte.x[n]

  • Processe os pontos para criar uma matriz de bits seguinte maneira. Aqui é um pequeno número escolhido para se adequar à sua noção de quão perto de uma linha reta você deseja pontos para cortar para. O critério será reconhecido pelos cognoscentos como exigindo que a linha reta através de e possua quase a mesma inclinação que a linha reta através de e .x[n]y[n]

    y[n]={1,if |(x[n+1]x[n])(x[n]x[n1])|<ϵ,0,otherwise.
    ϵ n + 1 ] )( n - 1 , x [ n - 1 ] ) ( n , x [ n ] ) ( n , x [ n ] ) ( n + 1 , x [x[n1],x[n],x[n+1](n1,x[n1])(n,x[n])(n,x[n])(n+1,x[n+1])
  • Se é uma matriz de dez ou corre tão alongados de s separados por corridas de s com vadios ocasional s aqui e ali para estragar a beleza, relaxe, você está no caminho certo. Caso contrário, se houver poucas execuções ou muitas execuções de s, repita a etapa anterior com um diferente .1 0 1 1 ϵy[n]1011ϵ

  • y[n]x[3]x[88]x[94]x[120]x[129], e assim por diante. Estenda A para a direita e B para a esquerda para descobrir onde eles se cruzam; estenda B para a direita e C para a esquerda para descobrir onde eles se cruzam etc. Parabéns, agora você tem um modelo linear contínuo e por partes para seus dados.

Dilip Sarwate
fonte
Roubou totalmente a minha resposta! =)
Phonon
Ideia interessante, mas, infelizmente, devido ao ruído no sinal, não consigo bons resultados.
P3trus
1
Essa expressão cujo magnituto está sendo comparado ao epsilon é na verdade uma aproximação à segunda derivada dos dados. Existem outras maneiras de calcular isso usando mais de três pontos que não respondem tanto ao ruído. Olhe Savitzky-Golay.
darenw
4

(Anos mais tarde) as funções lineares por partes são splines de grau 1, que podem ser solicitadas à maioria dos instaladores de spline. scipy.interpolate.UnivariateSpline, por exemplo, pode ser executado com k=1 um parâmetro de suavização s, com o qual você terá que brincar - consulte scipy-interpolation-with-univariate-splines .
No Matlab, veja como escolher nós .

Adicionado: encontrar nós ótimos não é fácil, porque pode haver muitos ótimos locais. Em vez disso, você atribui ao UnivariateSpline um destino s, soma do erro ^ 2, e permite que ele determine o número de nós. Após o ajuste, get_residual()obterá a soma real do erro ^ 2 e get_knots()os nós. Uma pequena mudança spode mudar bastante os nós, especialmente em ruídos altos - sim.
O gráfico mostra ajustes para uma função linear aleatória por partes + ruído para vários s.

Para ajustar constantes por partes, consulte Detecção de etapas . Isso pode ser usado para pw linear? Não sei; começar por diferenciar dados ruidosos aumentará o ruído errado.

Outras funções de teste e / ou links para documentos ou códigos seriam bem-vindos. Alguns links:
regressão linear por partes com nós como parâmetros
Splines lineares são muito sensíveis a onde os nós são colocados
seleção de nó para splines de regressão cúbica
Este é um problema complicado e a maioria das pessoas apenas seleciona os nós por tentativa e erro.
Uma abordagem que está crescendo em popularidade é usar splines de regressão penalizados.


Adicionado em março de 2014: a programação dinâmica é um método geral para problemas com subproblemas aninhados como este:

optimal k lines
    = optimal k - 1 lines up to some x
    + cost of the last line x to the end
over x  (all x in theory, nearby x in practice)

A programação dinâmica é muito inteligente, mas pode vencer a força bruta + heurística para esta tarefa?
Veja as excelentes notas de curso de Erik Demaine no MIT 6.006. Introdução aos algoritmos e
também à regressão linear segmentada pelo Google e à
síndrome de John Henry.


insira a descrição da imagem aqui

denis
fonte
O problema, pelo menos com scipy, é o posicionamento dos nós. O scipy usa nós igualmente espaçados.
P3trus
@ P3trus, sim, para começar, mas eles podem se mover - veja a trama. De qualquer forma, ele visa erro total, não nós.
Denis
@ P3trus Você já tentou usar o método de splines de regressão multivariada que seleciona automaticamente os pontos de interrupção iterativamente? cs.rtu.lv/jekabsons/regression.html
Atul Ingle
@Atul Ingle, um bom ponto de interrupção / seleção de nó é o mesmo problema, seja qual for o ajuste do spline. Se você conhece algum algoritmo diferente para o de pessoas de regressão / R, poderia postar um link, por favor?
Denis #
Está procurando pacotes no R / Matlab que fazem splines de regressão adaptável? Aqui: cran.r-project.org/web/packages/earth/index.html cran.r-project.org/web/packages/mda/index.html e também o ARESLab no Matlab para o qual eu já publiquei o link.
Atul Ingle
0

Pegue a derivada e procure por áreas de valor quase constante. Você precisaria criar o algoritmo para procurar por áreas com idealmente algum nível de +/- inclinação e isso daria a inclinação da linha para essa seção. Você pode querer realizar alguma suavização, como uma média deslizante, antes de fazer a classificação secional. O próximo passo seria obter a interseção y, que deve ser trivial nesse ponto.

Porten
fonte
derivado pode ser muito barulhento. Eu não acho que eu recomendaria isso.
Robert Bristow-johnson
0

Usar um filtro de tendência l1 é outra ideia:

Papel

Exemplo Online

SeanVN
fonte
1
Sua resposta é um pouco curta demais para ser construtiva! Por favor, considere fazer um esforço para expandi-lo de maneira pedagógica.
Sansuiso