传递熵

信息熵

由信息论之父——克劳德·艾尔伍德·香农提出，并首次用数学公式阐明了概率与信息冗余度之间的关系。为了方便理解，我们可以参照热力学里面的熵来理解信息熵。在热力学中，我们用熵来度量体系混乱程度，熵增代表物质向无规则方向发展即体系越来越混乱。其次，我们可以通过降温等手段来减熵即恢复有序性。同样，在信息论中，我们也可以这样理解：信息熵是用来度量体系不确定度的量，被观测的系统越不确定信息熵越大，系统越稳定信息熵越小。其次，我们可以通过被告知某个信息来降低系统的不确定性即‘减熵’。
但是值得注意的一点是，这个信息不一定可以降低不确定度，有时可能不会给你提供任何有效信息甚至可能会增大这个系统的不确定度。打个比方：玩投色子的游戏，点数1~6，投完色子让你去猜点数。在你猜之前呢，我告诉你的信息“这个点数是个位数”。那么对你来说这条信息就是无用的，他不会提高你猜对的概率。
下面是Shannon Entropy的公式：

我们可以看到信息熵是度量单个变量的不确定度的，我们在之后的运用就是建立在此基础上的。

因果关系

就是字面上的意思，获取原因的信息可以降低我们观测结果的不确定性。当我们知道了单个变量的信息熵之后，我们如何去度量变量之间是如何彼此影响的呢？我们应该如何去建立多个变量之间的因果关系呢？这个时候我们就需要引入另一个概念——传递熵(Tansfer Entropy)

传递熵

信息熵是用来研究变量与变量之间的信息的传递，可以计算这个信息传递能减少多少被观测系统的不确定度。当 X 对 Y 的传递熵 > Y 对 X 的传递熵时，我们就把 X 称为因，把 Y 称为果，并以此来建立两个变量之间的因果关系。
————————————————
版权声明：本文为优快云博主「Alphoseven」的原创文章，遵循CC 4.0 BY-SA版权协议，转载请附上原文出处链接及本声明。
原文链接：https://blog.youkuaiyun.com/Alphoseven/article/details/103491897

 

相信随着大数据的到来，会使得传递熵（transfer entropy）这种方法会受到更多人的关注，因为它是一种基于概率分布，香农熵，统计的方法得出时间序列间因果性的方法。这种方法首先提出是在PRL上由T.Schreiber 提出的，又用在生物系统中，所以这里大部分人可能不熟悉。由于transfer entropy 所需的时间序列长度较大，所以在普遍数据量较小的时代，只能用在神经信号和脑电图中。现在很多地方都意识到数据的重要性，各种传感器也被大量应用，原本不存在的数据来源慢慢的也被发掘了。相信transfer entropy 的春天就要来了。

什么是transfer entropy 它其实就是一个条件分布带来的探测到时间序列间的不对称性。说的学术一点：传递熵是在错误假设传递概率函数为p(in+1|in(k)),而不是p(in+1|in(k),jn(l))的情况下，预测系统状态额外需要的信息。这个信息由Y到X和由X到Y是不对称，这种不对称就带来了，驱动和响应的关系的建立。不过他和granger 因果性检验之间的等价性在一篇工作中已经证明。而且传递熵能对非线性时间序列应用，对这种granger的因果性也很敏感。

参见    http://blog.youkuaiyun.com/jjjphysics/article/details/24349289

作者：bioinfo2011
链接：https://www.jianshu.com/p/61b13fd053aa
来源：简书
著作权归作者所有。商业转载请联系作者获得授权，非商业转载请注明出处。

 

离散过程的传递熵表达式为：

y到x的传递熵实质为y的信息对于x不确定性大小的改变，即y传递给x的信息量的大小。因此，传递熵同样可以作为衡量因果性的指标。传递熵由于考虑的是变量间的信息量传递，而不需要假定变量间具有特定形式的关系，因此具有比Wiener-Granger因果性更好的适用性，尤其是对于非线性系统。

%X是一个时间序列，是一个一维列向量
%Y是一个时间序列，是一个一维列向量
%pieces是将空间离散多少分，一般在数据量足够大的情况下分的越多，能够得到的信息越多
%j是一个真实想要用的时间长度大小
function [en_1_2,en_2_1]=transfer_entropy(X,Y,pieces,j)
%% parpare some data that use in reconstruction phase space
d_x=[];    sit=size(X);
templ=randperm((sit(1,1)-1));
select=templ(1:j);
d_x(:,1)=X(select+1,1);d_x(:,2)=X(select,1);d_x(:,3)=Y(select+1,1);d_x(:,4)=Y(select,1);
X_max=max(X(:,1));X_min=min(X(:,1));Y_max=max(Y(:,1));Y_min=min(Y(:,1));
delta1=(X_max-X_min)/(2*pieces);delta2=(Y_max-Y_min)/(2*pieces);
%% calculate the data disturbution p(x(t)),p(x(t+1)),P(x(t+1),x(t)),p(x(t+1),x(t),y(t)),p(x(t),y(t)),p(y(t)),p(y(t+1)),p(y(t+1))
L1=linspace(X_min+delta1,X_max-delta1,pieces);L2=linspace(Y_min+delta2,Y_max-delta2,pieces);
dist1=zeros(pieces(1,1),2);
count=0;
for q1=1:pieces
    k1=L1(q1);k2=L2(q1);
    count=count+1;
    count1=0;count2=0;
    for i=1:j
        if d_x(i,2)>=(k1-delta1) && d_x(i,2)<=(k1+delta1)
            count1=count1+1;
        end
        if d_x(i,4)>=(k2-delta2) && d_x(i,4)<=(k2+delta2)
            count2=count2+1;
        end
    end
    dist1(count,1)=count1;dist1(count,2)=count2;
end

dist1(:,1)=dist1(:,1)/sum(dist1(:,1)); dist1(:,2)=dist1(:,2)/sum(dist1(:,2));

dist2=zeros(pieces(1,1),pieces(1,1),3);
for q1=1:pieces
    for q2=1:pieces
        k1=L1(q1);k2=L1(q2);
        k3=L2(q1);k4=L2(q2);
        count1=0;count2=0;count3=0;
        for i1=1:j
                if d_x(i1,1)>=(k1-delta1) && d_x(i1,1)<=(k1+delta1) && d_x(i1,2)>=(k2-delta1) && d_x(i1,2)<=(k2+delta1)
                    count1=count1+1;
                end
                if d_x(i1,3)>=(k3-delta2) && d_x(i1,3)<=(k3+delta2) && d_x(i1,4)>=(k4-delta2) && d_x(i1,4)<=(k4+delta2)
                    count2=count2+1;
                end
                if d_x(i1,2)>=(k1-delta1) && d_x(i1,2)<=(k1+delta1) && d_x(i1,4)>=(k4-delta2) && d_x(i1,4)<=(k4+delta2)
                    count3=count3+1;
                end
        end        
        dist2(q1,q2,1)=count1;dist2(q1,q2,2)=count2;dist2(q1,q2,3)=count3;        
    end
end
dist2(:,:,1)=dist2(:,:,1)/sum(sum(dist2(:,:,1)));dist2(:,:,2)=dist2(:,:,2)/sum(sum(dist2(:,:,2)));dist2(:,:,3)=dist2(:,:,3)/sum(sum(dist2(:,:,3)));

dist3=zeros(pieces(1,1),pieces(1,1),pieces(1,1),2);
for q1=1:pieces
    for q2=1:pieces
        for q3=1:pieces
            k1=L1(q1);k2=L1(q2);k3=L1(q3);
            k4=L2(q1);k5=L2(q2);k6=L2(q3);
            count1=0;count2=0;
            for i1=1:j
                if d_x(i1,1)>=(k1-delta1) && d_x(i1,1)<=(k1+delta1) && d_x(i1,2)>=(k2-delta1) && d_x(i1,2)<=(k2+delta1) && d_x(i1,4)>=(k6-delta2) && d_x(i1,4)<=(k6+delta2)  
                    count1=count1+1;
                end
                if d_x(i1,3)>=(k4-delta2) && d_x(i1,3)<=(k4+delta2) && d_x(i1,4)>=(k5-delta2) && d_x(i1,4)<=(k5+delta2) && d_x(i1,2)>=(k3-delta1) && d_x(i1,2)<=(k3+delta1)  
                    count2=count2+1;
                end
            end
            dist3(q1,q2,q3,1)=count1;dist3(q1,q2,q3,2)=count2;
        end
    end
end
dist3(:,:,:,1)=dist3(:,:,:,1)/sum(sum(sum(dist3(:,:,:,1))));dist3(:,:,:,2)=dist3(:,:,:,2)/sum(sum(sum(dist3(:,:,:,2))));

%% use the dosturbution calculate thansfer entropy
sum_f_1=0;sum_f_2=0;
for k1=1:pieces(1,1)
    for k2=1:pieces(1,1)
        if dist2(k1,k2,2)~=0 && dist1(k2,2)~=0  
           sum_f_1=sum_f_1-dist2(k1,k2,2)*log2(dist2(k1,k2,2)/dist1(k2,2));
        end
        
        if dist2(k1,k2,1)~=0 && dist1(k2,1)~=0
           sum_f_2=sum_f_2-dist2(k1,k2,1)*log2(dist2(k1,k2,1)/dist1(k2,1));
        end
    end
end
disp(sum_f_1);disp(sum_f_2)
sum_s_1=0;sum_s_2=0;
for k1=1:pieces(1,1)
    for k2=1:pieces(1,1)
        for k3=1:pieces(1,1)
            if dist3(k1,k2,k3,2)~=0 && dist2(k3,k2,3)~=0
               sum_s_1=sum_s_1-dist3(k1,k2,k3,2)*log2(dist3(k1,k2,k3,2)/dist2(k3,k2,3));
            end
            
            if dist3(k1,k2,k3,1)~=0 && dist2(k2,k3,3)~=0
               sum_s_2=sum_s_2-dist3(k1,k2,k3,1)*log2(dist3(k1,k2,k3,1)/dist2(k2,k3,3));
            end
        end
    end
end
disp(sum_s_1);disp(sum_s_2)

en_1_2=sum_f_1-sum_s_1;
en_2_1=sum_f_2-sum_s_2;
end
————————————————
版权声明：本文为优快云博主「jjjphysics」的原创文章，遵循CC 4.0 BY-SA版权协议，转载请附上原文出处链接及本声明。
原文链接：https://blog.youkuaiyun.com/jjjphysics/article/details/24349289