Fixed proper data supply to the KSG estimation for Schreiber heart-breath example

This commit is contained in:
joseph.lizier 2014-07-28 13:48:06 +00:00
parent 01291c883a
commit 2ac9c4b8e5
2 changed files with 15 additions and 9 deletions

View File

@ -0,0 +1,2 @@
See the octave folder for matlab demos (Octave is an open-source clone of Matlab).

View File

@ -6,7 +6,7 @@
% 22/1/2014
%
% Used to explore information transfer in the heart rate / breath rate example of Schreiber --
% but estimates TE using Kraskov-Grassberger estimation.
% but estimates TE using Kraskov-Stoegbauer-Grassberger estimation.
%
%
% Inputs
@ -42,7 +42,7 @@ function [teHeartToBreath, teBreathToHeart] = runHeartBreathRateKraskov(kHistory
fprintf('TE for heart rate <-> breath rate for Kraskov estimation with %d samples:\n', timeSteps);
% Using a single conditional mutual information calculator is the least biased way to run this:
% Using a KSG estimator for TE is the least biased way to run this:
teCalc=javaObject('infodynamics.measures.continuous.kraskov.TransferEntropyCalculatorKraskov');
for knnIndex = 1:length(knns)
@ -52,27 +52,31 @@ function [teHeartToBreath, teBreathToHeart] = runHeartBreathRateKraskov(kHistory
% Perform calculation for heart -> breath (lag 1)
teCalc.initialise(kHistory,1,lHistory,1,1);
teCalc.setProperty('k', sprintf('%d',knn));
teCalc.setObservations(octaveToJavaDoubleMatrix(heart), ...
octaveToJavaDoubleMatrix(chestVol));
teCalc.setObservations(octaveToJavaDoubleArray(heart), ...
octaveToJavaDoubleArray(chestVol));
teHeartToBreath(knnIndex) = teCalc.computeAverageLocalOfObservations();
% Perform calculation for breath -> heart (lag 1)
teCalc.initialise(kHistory,1,lHistory,1,1);
teCalc.setProperty('k', sprintf('%d',knn));
teCalc.setObservations(octaveToJavaDoubleMatrix(chestVol), ...
octaveToJavaDoubleMatrix(heart));
teCalc.setObservations(octaveToJavaDoubleArray(chestVol), ...
octaveToJavaDoubleArray(heart));
teBreathToHeart(knnIndex) = teCalc.computeAverageLocalOfObservations();
fprintf('TE(k=%d): heart->breath = %.3f, breath->heart = %.3f\n', knn, teHeartToBreath(knnIndex), teBreathToHeart(knnIndex));
fprintf('TE(k=%d): heart->breath = %.3f, breath->heart = %.3f nats\n', knn, teHeartToBreath(knnIndex), teBreathToHeart(knnIndex));
end
tElapsed = toc(starttime);
fprintf('Total runtime was %.1f sec\n', tElapsed);
plot(knns, teHeartToBreath, 'rx');
plot(knns, teHeartToBreath, 'rx', 'markersize', 10);
hold on;
plot(knns, teBreathToHeart, 'mo');
plot(knns, teBreathToHeart, 'mo', 'markersize', 10);
hold off;
legend(['TE(heart->breath)'; 'TE(breath->heart)']);
set (gca,'fontsize',26);
xlabel('K nearest neighbours', 'FontSize', 44, 'FontWeight', 'bold');
ylabel('TE', 'FontSize', 44, 'FontWeight', 'bold');
print('results.eps', '-deps', '-color');
end