ICA - Independência estatística e valores próprios da matriz de covariância

14

Atualmente, estou criando sinais diferentes usando o Matlab, misturando-os multiplicando-os por uma matriz de mixagem A e, em seguida, tentando recuperar os sinais originais usando o FastICA .

Até agora, os sinais recuperados são realmente ruins quando comparados aos originais, o que não era o que eu esperava.

Estou tentando ver se estou fazendo algo errado. Os sinais que estou gerando são os seguintes:

s1 = (-x.^2 + 100*x + 500) / 3000; % quadratic
s2 = exp(-x / 10); % -ve exponential
s3 = (sin(x)+ 1) * 0.5; % sine
s4 = 0.5 + 0.1 * randn(size(x, 2), 1); % gaussian
s5 = (sawtooth(x, 0.75)+ 1) * 0.5; % sawtooth

Sinais originais

Uma condição para o sucesso da ICA é que no máximo um sinal seja gaussiano, e observei isso na minha geração de sinal.

No entanto, outra condição é que todos os sinais sejam estatisticamente independentes.

Tudo o que sei é que isso significa que, dados dois sinais A e B, conhecer um sinal não fornece nenhuma informação em relação ao outro, ou seja: P (A | B) = P (A) onde P é a probabilidade .

Agora, minha pergunta é a seguinte: meus sinais são estatisticamente independentes? Existe alguma maneira de determinar isso? Talvez alguma propriedade que deva ser observada?

Outra coisa que notei é que, quando calculo os valores próprios da matriz de covariância (calculados para a matriz que contém os sinais mistos), o espectro eigens parece mostrar que existe apenas um componente principal (principal) . O que isso realmente significa? Não deveria haver 5, já que tenho 5 (supostamente) sinais independentes?

Por exemplo, ao usar a seguinte matriz de mistura:

A =

0.2000    0.4267    0.2133    0.1067    0.0533
0.2909    0.2000    0.2909    0.1455    0.0727
0.1333    0.2667    0.2000    0.2667    0.1333
0.0727    0.1455    0.2909    0.2000    0.2909
0.0533    0.1067    0.2133    0.4267    0.2000

Os autovalores são: 0.0000 0.0005 0.0022 0.0042 0.0345(apenas 4!)

Quando se utiliza a matriz identidade como a matriz de mistura (isto é, os sinais mistos são as mesmas que as originais), o eigenspectrum é: 0.0103 0.0199 0.0330 0.0811 0.1762. Ainda existe um valor muito maior que o resto.

Obrigado pela ajuda.

Peço desculpas se as respostas às minhas perguntas são dolorosamente óbvias, mas sou realmente novo em estatísticas, ICA e Matlab. Obrigado novamente.

EDITAR

Eu tenho 500 amostras de cada sinal, no intervalo [0,2, 100], nas etapas de 0,2, ou seja, x = 0: 0,1: 100.

Além disso, dado o modelo ICA: X = As + n (não estou adicionando nenhum ruído no momento), estou me referindo ao espectro eigens da transposição de X, ou seja, eig (cov (X ')).

ATUALIZAR

Conforme sugerido (consulte os comentários), tentei o FastICA em apenas 2 sinais. Os resultados foram muito bons (veja a foto abaixo). A matriz de mistura utilizada foi A = [0.75 0.25; 0.25 0.75]. No entanto, o espectro eigens 0.1657 0.7732ainda mostrava apenas um componente principal principal.

Minha pergunta, portanto, se resume ao seguinte: Que função / equação / propriedade posso usar para verificar se vários vetores de sinais são estatisticamente independentes?

Seno e Gaussiano - FastICA

Rachel
fonte
1
Excelente pergunta. Eu perguntei sobre como podemos saber quando dois sinais são independentes aqui ( dsp.stackexchange.com/questions/1242/… ), mas não fui muito longe com isso. :-) Eu também sou novo na ICA, mas talvez eu possa esclarecer um pouco.
Spacey
@Mohammad Você ainda está interessado em uma resposta a essa pergunta? De bom grado colocarei uma recompensa nele para atrair interesse.
Phonon
@Mohammad eu votei na sua pergunta. Espero que você tenha uma boa resposta, está realmente relacionada à minha. Eu tenho lido os comentários até agora, e há muitas estatísticas acontecendo que eu não entendo. Você conseguiu encontrar uma maneira definitiva de concluir se dois sinais são independentes ou não?
Rachel
@ Rachel Ainda não no momento, mas vou pesquisar um pouco mais e informá-lo. É um conceito muito importante, que eu sinto que geralmente é vitrificado, infelizmente.
Spacey
Obrigado @Mohammad. Concordo. Sinais independentes observam a propriedade E (s1, s2) = E (s1) x E (s2), mas não sei como calculá-lo para sinais reais.
22412 Rachel

Respostas:

8

Os sinais 3 e 5 parecem estar bastante correlacionados - eles compartilham seu primeiro harmônico. Se eu recebesse duas misturas dessas, não seria capaz de separá-las, ficaria tentado a colocar a harmônica comum como um sinal e as harmônicas mais altas como segundo sinal. E eu estaria errado! Isso pode explicar o valor próprio ausente.

Os sinais 1 e 2 também não parecem independentes.

Uma "verificação de sanidade" rápida e suja para a independência de duas séries é fazer um gráfico (x, y) de um sinal contra o outro:

plot (sig3, sig5)

e depois fazer o mesmo gráfico (x, y) com um sinal embaralhado:

indices = randperm(length(sig3))
plot(sig3(indices), sig5)

Se os dois gráficos tiverem uma aparência diferente, seus sinais não serão independentes. De maneira mais geral, se o gráfico (x, y) dos dados mostra "características", desimetrias etc., é um mau presságio.

Testes adequados de independência (e essas são as funções objetivas usadas no loop de otimização da ACI) incluem, por exemplo, informações mútuas.

A ICA está recuperando os sinais mais independentes, cuja mistura linear produz seus dados de entrada . Ele funcionará como um método de separação de sinal e recuperará os sinais originais apenas se eles fossem maximamente independentes, de acordo com o critério de otimização usado na sua implementação de ICA.

pichenettes
fonte
1
Pergunta: Se os 5 sinais no caso dela fossem de fato todos independentes, esperaríamos que NÃO houvesse componentes principais corretos? (Em outras palavras, todos os valores próprios seriam iguais). Geometricamente, teríamos uma 'nuvem' guassiana em 5 dimensões, concorda?
Spacey
Eu também havia contatado um autor da ACI sobre a remoção de dois sinusóides de uma mistura, e ele disse que isso poderia ser feito com a ACI. Isso me confunde um pouco com base no que você está dizendo em relação aos sinais 3 e 5 porque (eu concordo com você) eles parecem correlacionados.
Spacey
@ pichenettes Eu plotei esses gráficos como você sugeriu - e os gráficos realmente têm uma aparência diferente. Infelizmente ainda estou preso a como testar a independência. Eu realmente preciso de uma maneira de gerar sinais estatisticamente independentes, para que eu possa avaliar o desempenho do FastICA.
22412 Rachel
x1[n]x2[n]
@Mohammad Eu não gravei minha própria voz, mas tentei usar o FastICA em uma mistura de sinais sinusodiais e gaussianos. Eu acho que eles são independentes .. FastICA teve um bom desempenho, mas o espectro de eigens ainda era estranho. Vou atualizar minha pergunta para mostrar os resultados.
244 Rachel
7

Não sou especialista em ACI, mas posso falar um pouco sobre independência.

Como alguns dos comentários mencionaram, a independência estatística entre duas variáveis ​​aleatórias pode ser mais ou menos interpretada como "a quantidade de informação que a observação de uma variável fornece sobre outra".

XYXYp(x,y)XYp(x,y)=p(x)p(y)

p(x,y)

XYXYp(X=i,Y=j)=pijP(X=i)=piP(Y=j)=pj

Eu(X,Y)=EujpEujregistropEujpEupj

Aqui está um código matlab que irá gerar dois sinais independentes de uma distribuição conjunta construída e dois de uma distribuição conjunta não independente e depois calcular as informações mútuas das juntas.

A função "computeMIplugin.m" é uma função simples que escrevi que calcula as informações mútuas usando a fórmula de soma acima.

Ndist = 25;
xx = linspace(-pi, pi, Ndist);

P1 = abs(sin(xx)); P2 = abs(cos(xx)); 
P1 = P1/sum(P1); P2 = P2/sum(P2); % generate marginal distributions

%% Draw independent samples.
Nsamp = 1e4;
X = randsample(xx, Nsamp, 'true', P1);
Y = randsample(xx, Nsamp, 'true', P2);

Pj1 = P1'*P2;
computeMIplugin(Pj1)

% I get approx 8e-15 ... independent!

% Now Sample the joint distribution 
cnt = {}; cnt{1} = xx; cnt{2} = xx; % bin centers
Pj1_samp= hist3([X' Y'],cnt); Pj1_samp = Pj1_samp/sum(Pj1_samp(:));
computeMIplugin(Pj1_samp)
% I get approx .02; since we've estimated the distribution from
% samples, we don't know the true value of the MI. This is where
% a confidence interval would come in handy. We'd like to know 
% whether value of MI is significantly different from 0. 

% mean square difference between true and sampled?
% (this is small for these parameter settings... 
% depends on the sample size and # bins in the distribution).
mean( (Pj1_samp(:) - Pj1(:)).^2)

%% Draw samples that aren't independent. 

tx = linspace(0,30,Nsamp);
X = pi*sin(tx);
Y = pi*cos(tx);

% estimate the joint distribution
cnt = {}; cnt{1} = xx; cnt{2} = xx; % bin centers
Pj2= hist3([X' Y'],cnt); Pj2 = Pj2/sum(Pj2(:));
computeMIplugin(Pj2)

% I get 1.9281  - not independent!

%% make figure
figure(1); 
colormap gray
subplot(221)
imagesc(xx,xx,Pj1_samp)
title('sampled joint distribution 1')
subplot(222)
imagesc(xx,xx,Pj2)
title('sampled joint distribution 2')
subplot(223)
imagesc(xx,xx,Pj1)
title('true joint distribution 1')

Mais uma vez, isso pressupõe que você tenha uma boa estimativa da distribuição conjunta (junto com outras suposições mais alegres), mas deve ser útil como regra geral.

Sim
fonte
Essa é uma boa resposta sydeulissie obrigado, vou ter que olhar um pouco mais fundo.
Spacey
Primeiro de tudo, obrigado pela resposta longa, foi muito informativo. Eu só tenho algumas perguntas. Você mencionou o uso de um teste qui-quadrado. Eu olhei para ele e parece realmente interessante, mas como posso usá-lo em sinais? Não pode ser aplicado apenas a dados categóricos?
Rachel
Além disso, você está usando Pj1 = P1 '* P2 para calcular a distribuição da junta, correto? Mas, tecnicamente, acredito que isso não pode ser feito. Talvez você esteja fazendo isso porque está assumindo que os sinais originais são independentes e o resultado é válido? Mas como você pode calcular as informações mútuas - já que o resultado depende da distribuição conjunta ..? Pode ser que eu tenha entendido algo errado, mas gostaria de um esclarecimento, por favor.
Rachel
Ficarei feliz em - embora demore um pouco antes de eu ter tempo :).
sim
Obrigado @sydeulissie. Gostaria de uma resposta que não assuma que tenho conhecimento da distribuição conjunta, por favor.
Rachel
3

Como mencionado acima, ambos os sinais 3 e 5 parecem estar bastante correlacionados e têm um período semelhante.

Podemos pensar em dois sinais sendo correlacionados se pudermos mudar uma das fontes para a esquerda ou direita e aumentar ou diminuir sua amplitude para que ela se encaixe em cima da outra fonte. Observe que não estamos alterando a frequência da fonte, apenas realizando uma mudança de fase e amplitude.

No caso acima, podemos mudar a fonte 3 para que seus picos coincidam com a fonte 5. Esse é o tipo de coisa que interferirá na extração da fonte ao usar o ICA devido à suposição de independência.

Nota : Uma boa ilustração do conceito acima é pensar em duas ondas sinusoidais. Ambos são completamente determinísticos. Se ambos tiverem a mesma frequência (mesmo com fases diferentes), estarão perfeitamente correlacionados e a ACI não poderá separá-los. Se, em vez disso, eles tiverem frequências diferentes (que não são múltiplos inteiros um do outro), eles serão independentes e poderão ser separados.

Abaixo está um código Matlab para você ver por si mesmo

%Sine waves of equal frequency
X = 1:1000;
Y(1,:) = sin(2*pi*X*10/1000);
Y(2,:) = sin(1+2*pi*X*10/1000);

figure
subplot(3,2,1)
plot(Y(1,:))
title('Initial Source 1')
subplot(3,2,2)
plot(Y(2,:))
title('Initial Source 2')
A = [1, 2; 4, -1];
Y = A*Y;
subplot(3,2,3)
plot(Y(1,:))
title('Signal 1')
subplot(3,2,4)
plot(Y(2,:))
title('Signal 2')

Z = fastica(Y);

subplot(3,2,5)
plot(Z(1,:))
title('Source 1')
subplot(3,2,6)
plot(Z(2,:))
title('Source 2')

%Sine waves of different frequency
X = 1:1000;
Y(1,:) = sin(2*pi*X*10/1000);
Y(2,:) = sin(1+2*pi*X*8/1000);

figure
subplot(3,2,1)
plot(Y(1,:))
title('Initial Source 1')
subplot(3,2,2)
plot(Y(2,:))
title('Initial Source 2')
A = [1, 2; 4, -1];
Y = A*Y;
subplot(3,2,3)
plot(Y(1,:))
title('Signal 1')
subplot(3,2,4)
plot(Y(2,:))
title('Signal 2')

Z = fastica(Y);

subplot(3,2,5)
plot(Z(1,:))
title('Source 1')
subplot(3,2,6)
plot(Z(2,:))
title('Source 2')

Observe que, para as ondas da mesma frequência, o ICA apenas retorna os sinais de entrada, mas para diferentes frequências, ele retorna as fontes originais.

rwolst
fonte
2

Rachel,

Até agora, em minhas pesquisas, consegui encontrar algo chamado ' Teste Qui-Quadrado de Independência ', mas não sei ao certo como funciona no momento, mas pode valer a pena dar uma olhada.

Spacey
fonte
Encontrei esses dois tutoriais explicando como executar o teste do qui-quadrado: ling.upenn.edu/~clight/chisquared.htm & math.hws.edu/javamath/ryan/ChiSquare.html . No entanto, o teste pode ser realizado apenas com dados categóricos. Eu não sei se isso pode ser aplicado a nossas observações de sinal ..
Rachel