2012-08-22 4 views
1

MATLAB의 지수 커널 (k)과 함께 두 개의 스파이크 (Spike라고 함)가 포함 된 시계열을 컨볼루션하고 싶습니다. convolved response "calcium1"을 호출하십시오. 커널과 디콘 볼 루션 (deconvolution)을 사용하여 원래의 스파이크 ("reconSpike") 데이터를 복구하려고합니다. 다음 코드를 사용하고 있습니다.MATLAB의 deconv() 함수는 conv() 함수를 반전하지 않습니다.

k1=zeros(1,5000); 
k1(1:1000)=(1.1.^((1:1000)/100)-(1.1^0.01))/((1.1^10)-1.1^0.01); 
k1(1001:5000)=exp(-((1001:5000)-1001)/1000); 
k1(1)=k1(2); 

spike = zeros(100000,1); 
spike(1000)=1; 
spike(1100)=1; 

calcium1=conv(k1, spike); 
reconSpike1=deconv(calcium1, k1); 

문제 reconSpike의 끝에서, I 원본 데이터 아니었다 매우 큰 진폭 높은 파도의 덩어리를 얻을 수 있다는 것이다. 누구나 왜 고쳐야하는지 알고 있습니까?

감사합니다.

답변

1

스파이크 벡터를 k1 벡터와 동일한 길이로 유지하면 저에게 효과적입니다. 예 :

k1=zeros(1,5000); 
    k1(1:1000)=(1.1.^((1:1000)/100)-(1.1^0.01))/((1.1^10)-1.1^0.01); 
    k1(1001:5000)=exp(-((1001:5000)-1001)/1000); 
    k1(1)=k1(2); 

    spike = zeros(5000, 1); 
    spike(1000)=1; 
    spike(1100)=1; 

    calcium1=conv(k1, spike); 
    reconSpike1=deconv(calcium1, k1); 

당신이 다른 이유는 무엇입니까?

+0

나는 자유가 시간 [0 : 5000] 밖에서 스파이크를 일으키기를 원하기 때문에 나는 그것을 다르게 만들었다. 예를 들어, t = 6000에서 스파이크를 원한다면 그러나 5000보다 큰 모든 값이 0이기 때문에 나는 왜이 문제가 있는지 확신하지 못합니다. 아마도 k1의 길이를 늘리고 모두 0으로 설정하면 해결할 수 있습니다. – jfeng92

+0

예 그걸로 갈 것입니다. – Dan

+0

불행하게도, k1의 길이를 100000으로 늘리는 것이 효과가 없으며 여전히 이러한 파도가 발생합니다. 아마도 이것은 deconv 함수에 내재 된 결함입니다. 어떤 대안을 생각해 볼 수 있습니까? – jfeng92

1

MATLAB의 deconvolution 알고리즘 또는 부동 소수점 정밀 문제 (또는 둘 다)에 문제가 있습니다. 디콘 볼 루션 중에 발생하는 모든 디비전 및 뺄셈으로 인해 부동 소수점 정밀도가 의심 스럽지만 MathWorks에 직접 문의하여 가치를 묻는 것이 좋습니다. MATLAB 문서 당

, [q,r] = deconv(v,u) 있다면, v = conv(u,q)+r도 (즉, deconv의 출력을 항상 만족한다) 길게한다. 귀하의 경우 이것은 심각하게 위반됩니다. 스크립트의 끝에 다음을 넣어 :

[reconSpike1 rem]=deconv(calcium1, k1); 
max(conv(k1, reconSpike1) + rem - calcium1) 

내가 6.75e227를 얻을하지 제로 ;-) 다음 spike 6000의 길이를 변경하려고 할 것입니다; 작은 숫자 (~ 1e-15)를 얻을 것입니다. 점차적으로 길이를 늘리십시오 spike; 오류는 더 커지고 커질 것입니다. 스파이크에 0이 아닌 요소를 하나만 넣으면이 동작은 발생하지 않습니다. 오류는 항상 0입니다. 말된다; 모든 MATLAB은 모든 것을 같은 수로 나눕니다. 그것은 더 조언을 제공하기 어렵다, 그래서

maximum error for length(v) = 100 is 0.000000 
maximum error for length(v) = 1000 is 14.910770 

내가 당신이 정말로 달성하려는 모르겠어요 :

v = random('uniform', 1,2,100,1); 
u = random('uniform', 1,2,100,1); 
[q r] = deconv(v,u); 
fprintf('maximum error for length(v) = 100 is %f\n', max(conv(u, q) + r - v)) 
v = random('uniform', 1,2,1000,1); 
[q r] = deconv(v,u); 
fprintf('maximum error for length(v) = 1000 is %f\n', max(conv(u, q) + r - v)) 

출력은 다음과 같습니다

임의의 벡터를 사용하여 간단한 데모입니다 . 그러나 펄스가 쌓여서 각 펄스에 대한 정보를 추출하는 데 문제가있는 경우 이는 까다로운 문제 일 수 있습니다. 나는 이와 같은 일을하는 사람들을 알고 있습니다. 그래서 어떤 참고 문헌이 나에게 알려 지길 원한다면 나는 그들에게 물어볼 것입니다.

1

당신은 디컨 볼 루션은 단순히 컨볼 루션을 취소 할 수 있습니다 기대하지한다. 이것은 디콘 볼 루션이 ill-posed problem이기 때문입니다.

문제는 컨볼 루션이 정수 연산자 (연속적인 경우 사용자가 int f(x) g(x-t) dx 또는 유사한 것을 적어 두는 경우)에서 비롯된다. 이제, 적분 계산의 역 (- 컨볼 루션)은 차별화를 적용하는 것입니다. 불행히도, 차동 입력의 잡음을 증폭시킵니다. 따라서 통합에 약간의 오류 만있는 경우 (부동 소수점 부정확성으로 이미 충분할 수도 있음) 차별화 된 결과가 달라질 수 있습니다.

이러한 증폭을 완화 할 수있는 방법이 있지만 응용 프로그램별로 시도해야합니다.

관련 문제