Suavização de dados 2D

8

Os dados consistem em espectros ópticos (intensidade da luz em relação à frequência) gravados em momentos variados. Os pontos foram adquiridos em uma grade regular em x (tempo), y (frequência). Para analisar a evolução do tempo em frequências específicas (um aumento rápido, seguido de uma deterioração exponencial), eu gostaria de remover parte do ruído presente nos dados. Esse ruído, para uma frequência fixa, provavelmente pode ser modelado como aleatório com distribuição gaussiana. Em um horário fixo, no entanto, os dados mostram um tipo diferente de ruído, com picos espúrios grandes e oscilações rápidas (+ ruído gaussiano aleatório). Tanto quanto posso imaginar, o ruído ao longo dos dois eixos não deve ser correlacionado, pois tem origens físicas diferentes.

Qual seria um procedimento razoável para suavizar os dados? O objetivo não é distorcer os dados, mas remover artefatos ruidosos "óbvios". (e a suavização excessiva pode ser ajustada / quantificada?) Não sei se a suavização em uma direção independentemente da outra faz sentido ou se é melhor suavizar em 2D.

Eu li coisas sobre estimativa de densidade de kernel 2D, interpolação de polinômio / spline 2D, etc., mas não estou familiarizado com o jargão ou a teoria estatística subjacente.

Eu uso R, para o qual vejo muitos pacotes que parecem relacionados (MASS (kde2), campos (smooth.2d), etc.), mas não consigo encontrar muitos conselhos sobre qual técnica aplicar aqui.

Fico feliz em saber mais, se você tiver referências específicas para me indicar (ouvi dizer que MASS seria um bom livro, mas talvez seja muito técnico para um não estatístico).

Edit: Aqui está um espectrograma fictício representativo dos dados, com fatias ao longo das dimensões de tempo e comprimento de onda.

image2d

O objetivo prático aqui é avaliar a taxa de decaimento exponencial no tempo para cada comprimento de onda (ou caixas, se muito barulhento).

baptiste
fonte
Quantas frequências foram realizadas as medições? Se não for um número grande, seria prático ver isso como um conjunto de séries temporais individuais (mas relacionadas), uma para cada frequência específica?
Peter Ellis
@PeterEllis um grande número (digamos 500, mas por uma questão de generalidade poderia ser ainda maior)
Baptiste
Meu palpite é tratá-los como mais de 500 séries temporais correlacionadas e usar técnicas de séries temporais como média móvel ou suavização exponencial; e use somente a suavização 2d depois e somente se necessário para uma representação gráfica estilizada. No entanto, não tenho o suficiente para fazer disso uma "resposta" adequada.
Peter Ellis
1
Eu procuraria métodos 'robustos'. Esses métodos tentam reduzir o peso de valores discrepantes. Por exemplo, existe um algoritmo de spline robusto em R.
user12719
existe algo específico para a variável temporal que faz da série temporal um tipo específico de análise estatística?
precisa

Respostas:

4

Você precisa especificar um modelo que separa o sinal do ruído.

Existe o componente de ruído no nível de medição que você assume gaussiano. Os outros componentes, dependentes nas medições:

  • "Esse ruído, para uma frequência fixa, provavelmente pode ser modelado como aleatório com distribuição gaussiana". Precisa de esclarecimentos - o componente de ruído é comum a todos os pontos do tempo, dada a frequência? O desvio padrão é o mesmo para todas as frequências? Etc.

  • "Em um horário fixo, no entanto, os dados mostram um tipo diferente de ruído, com grandes picos espúrios e oscilações rápidas" Como você separa isso do sinal, pois provavelmente você está interessado na variação da intensidade através da frequência. A variação interessante é de alguma forma diferente da variação desinteressante? Em caso afirmativo, como?

Oscilações espúrias ou ruído não gaussiano em geral não são um grande problema, se você tiver uma idéia realista de suas características. Pode ser modelado transformando os dados (e depois usando um modelo gaussiano) ou explicitamente usando uma distribuição de erro não gaussiana. O ruído de modelagem correlacionado sobre as medições é mais desafiador.

Dependendo de qual é o seu modelo de ruído e dados, você pode modelar os dados com uma ferramenta de uso geral, como os GAMs no pacote mgcv, ou pode precisar de uma ferramenta mais flexível, que leve facilmente a uma configuração bayesiana bastante personalizada . Existem ferramentas para esses modelos, mas se você não for estatístico, aprender a usá-los levará um tempo.

Eu acho que ou uma solução específica para análise espectral ou o pacote mgcv são suas melhores apostas.

Scellus
fonte
Bom conselho, obrigado, preciso dar uma olhada nessas opções e pensar com mais cuidado na descrição do ruído.
baptiste
1
O ruído nos espectros ópticos geralmente depende da intensidade da luz medida ("contar" os fótons é um processo de Poisson) e geralmente também do comprimento de onda / frequência (devido às características do detector). Há muitos processos contribuindo com ruído instrumental, veja, por exemplo, Skoog & Leary: Principles of Instrumental Analysis. O tipo predominante de ruído dependerá do tipo de instrumento (e do experimento). O gráfico d ao longo do tempo mostra uma clara dependência da magnitude, sugerindo que Baptiste tem medições de intensidade (em oposição a, por exemplo, espectros de absorvância).
cbeleites infeliz com SX
2

Uma série temporal de espectros me sugere um experimento de cinética , e há uma quantidade bem estabelecida de literatura quimiométrica sobre isso.

O que você sabe sobre os espectros? Que tipo de espectro são eles? Você pode razoavelmente esperar que você tenha apenas duas espécies, eduto e produto?

Você pode razoavelmente assumir a bilinearidade, ou seja, os espectros medidos em um determinado momento são uma combinação linear das concentrações do componente vezes o espectro puro do componente : XCS

X(nspc×nWeu)=C(nspc×ncomp)S(ncomp×nWeu)

Você diz que deseja estimar uma deterioração exponencial (nas concentrações). Isso, juntamente com a bilinearidade, sugere-me a resolução de curvas multivariadas (MCR). Essa é uma técnica que permite que você use as informações que você possui (por exemplo, espectros de componentes puros de algumas substâncias ou suposições sobre o comportamento da concentração como decaimento exponencial) durante o ajuste do modelo.

Até onde eu sei, é bastante comum suavizar as concentrações de acordo com algum modelo, por exemplo, cinético, mas é muito menos comum suavizar os espectros. No entanto, o algoritmo permite fazer isso. No verão, perguntei a Anna se elas impõem restrições de suavidade, mas ela me disse que não (e os bons espectroscopistas odeiam a suavização em vez de medir bons espectros ;-)). Muitas vezes, isso também não é necessário, porque a agregação de informações de todos os espectros já produzirá boas estimativas dos espectros de componentes puros.

Eu suavizei "espectros de componentes" (de fato, componentes principais) duas vezes recentemente ( Dochow et al .: Dispositivo Raman no chip e fibras de detecção com grade de Bragg de fibra para análise de soluções e partículas, LabChip, 2013 e Dochow el al. : Chip microfluídico de quartzo para identificação de células tumorais por espectroscopia Raman em combinação com armadilhas ópticas, AnalBioanalChem, aceito), mas nesses casos meu conhecimento espectroscópico me dizia que eu posso fazer isso. Aplico regularmente uma interpolação de downsampling e suavização em meus espectros Raman ( hyperSpec::spc.loess).

Como saber o que é muita suavização? Eu acho que a única resposta possível é "conhecimento especializado sobre o tipo de espectroscopia e experimento".


editar: reli a pergunta e você diz que deseja estimar a deterioração em cada comprimento de onda. No entanto, isso é verdade ou você deseja estimar a deterioração de diferentes espécies com espectros sobrepostos?

cbeleites descontentes com o SX
fonte
Obrigado pelas referências. Embora a amostra realmente não tenha duas espécies, é um pouco semelhante (dois processos físicos distintos para distinguir). Vou dar uma olhada mais de perto quando voltar de uma conferência.
baptiste
@ batiste: Tenha uma boa conferência. Você se importa em dizer que tipo de processos você tem? Ou seja, você pode assumir que "dentro" de cada processo as características espectrais são as mesmas ou as oscilações podem "se mover" sobre o espectro (a frequência é ambígua se você tiver um padrão de oscilação no espectro )?
cbeleites infeliz com SX
1

Os dados consistem em espectros ópticos (intensidade da luz em relação à frequência) gravados em> tempos variados. Os pontos foram adquiridos em uma grade regular em x (tempo), y (frequência).

EuntensEuty=f(tEume,freqvocêency)fsendo uma soma de funções básicas (por exemplo, splines-b) e coeficientes. Um conjunto limitado de funções básicas reduz diretamente a rugosidade e, portanto, cancela boa parte do ruído branco.

Eu li coisas sobre estimativa de densidade de kernel 2D, interpolação de polinômio / spline 2D, etc.

...

Eu uso R, para o qual vejo muitos pacotes que parecem relacionados (MASS (kde2), campos (smooth.2d), etc.), mas não consigo encontrar muitos conselhos sobre qual técnica aplicar aqui.

Você mencionou a interpolação de spline, mas não mencionou o pacote fda, que implementa de maneira bastante fácil e fácil a expansão da função básica que mencionei acima. O conjunto de medições simultâneas de tempo, frequência e intensidade (ordenadas como uma matriz tridimensional) pode ser capturado como um objeto de dados funcionais bivariados, ver. por exemplo, a função 'Data2fd'. Além disso, vários procedimentos de suavização estão disponíveis na embalagem, todos projetados para cancelar o ruído branco ou a "rugosidade" nas medições de processos inerentemente suaves.

O artigo da Wikipedia descreve o problema do ruído branco no FDA da seguinte maneira:

Os dados podem ser tão precisos que o erro pode ser ignorado, pode estar sujeito a um erro de medição substancial ou até ter um relacionamento indireto complexo com a curva que eles definem. ... registros diários de precipitação em uma estação meteorológica são tão variáveis ​​que exigem análises cuidadosas e sofisticadas para extrair algo como uma curva de precipitação média.

O FDA fornece as ferramentas para esses casos. Isso não se traduz no seu caso?

... mas não estou familiarizado com o jargão ou a teoria estatística subjacente ...

... mas não consigo encontrar muitos conselhos sobre qual técnica aplicar aqui ...

Com relação à fda: eu não era nenhum dos dois, mas o livro de Ramsay e Silverman sobre a FDA (2005) torna o básico muito bem acessível e Ramsay Hooker e Graves (2009) traduzem diretamente os insights do livro em código R. Ambos os volumes devem estar disponíveis como e-books em uma biblioteca universitária de estatística, biociências, climatologia ou psicologia. O Google também exibirá mais alguns links que não posso postar aqui.

Lamento não poder fornecer uma solução mais direta para o seu problema. No entanto, o FDA me ajudou muito depois que eu descobri para que serve.

user1966337
fonte
Isso é útil, obrigado. Eu esperava ouvir uma perspectiva mais global do que apenas uma técnica em particular, mas se for essa a que devo usar, é tudo de bom.
precisa
Obrigado pelo crédito. Por fim, ninguém, a não ser você ou seus colegas imediatos, poderá decidir qual é a metodologia apropriada. Mas, à luz do que você descreveu, eu daria uma olhada no FDA em geral. Isso pode lhe dar mais algumas idéias de como analisar seus dados.
user1966337
@ user1966337: FYI, em espectroscopia óptica, as intensidades em comprimentos de onda distintos geralmente têm um significado distinto, para que você possa tratá-las como variáveis ​​de um modelo bilinear (fisicamente significativo) com poucos componentes, levando a um modelo mais restrito dos dados. Às vezes, porém, você tem efeitos que não permitem isso e onde o FDA seria mais apropriado.
cbeleites descontente com SX
1

Sendo um físico simples, não um especialista em estatística, eu adotaria uma abordagem simples. As duas dimensões são de naturezas diferentes. Faria sentido suavizar ao longo do tempo com um algoritmo e suavizar ao longo do comprimento de onda com outro.

Os algoritmos reais que eu usaria: para comprimento de onda, Savitzky-Golay com uma ordem superior, 6 talvez 8.

Com o tempo, se esse exemplo for típico, esse salto repentino e um declínio mais ou menos exponencial o tornam complicado. Eu tive dados experimentais e imagens barulhentas, exatamente assim. Se métodos simples e diretos não ajudarem o suficiente, tente um Gaussiano mais suave, mas suprima seu efeito próximo ao salto, como detectado por um detector de borda. Suavize e amplie a saída do detector de bordas, normalize-a para ir de 0,0 a 1,0 e use-a para selecionar entre a imagem original e a imagem suavizada de Gauss, pixel por pixel.

DarenW
fonte
0

@ batiste: Fico feliz que você adicionou o enredo como sugeri. Isso ajuda muito:

Portanto, se eu entendi direito, seu objetivo prático é avaliar a taxa de decaimento exponencial para cada comprimento de onda; então vamos fazer exatamente isso! Defina uma função que você deseja minimizar para cada comprimento de onda separadamente e minimize-a.

Vejamos um único comprimento de onda, como no gráfico inferior direito.

τ

τ^=umargmEunτtEu||e-tEu/τ-dEu||2

ττ

Posteriormente, se você acredita que o comprimento de onda adjacente deve ter constantes de decaimento semelhantes, você pode incorporar isso a um critério de otimização mais elaborado.

Se alguma coisa, eu sugiro que você leia um livro de otimização que deve ser lido: a otimização convexa de Boyd .

Espero que isto ajude!

zorbar
fonte
desculpe, mas parece haver um mal-entendido: eu estou familiarizado com a otimização não-linear; aqui, eu gostaria de saber quais técnicas de suavização eu posso usar nesses dados quando o ajuste em cada comprimento de onda não é confiável por causa do ruído nas duas dimensões. É verdade que meu exemplo idiota parece bastante viável, mas se eu tivesse acrescentado mais ruído, teria sido mais difícil de visualizar. Gosto da abordagem do FDA sugerida anteriormente, pois abrange tanto a parte do encaixe quanto a suavização em uma metodologia.
precisa saber é