copent包0.2版已在CRAN发布

是 否 +2 论坛币

k人 参与回答

经管之家送您一份

应届毕业生专属福利!

求职就业群

赵安豆老师微信：zhaoandou666

经管之家联合CDA

送您一个全额奖学金名额~ !

立即领取

感谢您参与论坛问题回答

经管之家送您两个论坛币！

+2 论坛币

copent包0.2版已在CRAN发布，同时在github发布。可从CRAN或github安装。
原来的0.1版提供copula熵估计，可进行统计独立性测试。
0.2版又增加了基于copula熵估计的条件独立性测试和传递熵（Transfer Entropy）估计两个功能，分别由 ci 和 transent 两个函数实现。
欢迎大家使用，并提出宝贵意见。

CRAN:
https://cran.r-project.org/web/packages/copent/index.html
GitHub:
https://github.com/majianthu/copent

此包的python版本同时更新到0.2版，增加了同样功能的同名函数。已在PyPI上发布。

扫码加我 拉你入群

请注明：姓名-公司-职位

以便审核进群资格，未注明则拒绝

分享0 收藏0 回帖

关键词：Open COPE CRAN CRA Pen

一定支持， 我可以下载源代码看不？
通过 git for windows ，
用 Rstudio 版本控制下载项目；

check , build  ， 看看代码。

源代码是这个。

1. ##################################################################################
   
2. ###  Estimating Copula Entropy and Transfer Entropy
   
3. ###  2021-03-20
   
4. ###  by MA Jian (Email: majian03@gmail.com)
   
5. ###
   
6. ###  Parameters
   
7. ###   x                : N * d data, N samples, d dimensions
   
8. ###   k         : kth nearest neighbour, parameter for kNN entropy estimation
   
9. ###   dt        : distance type [1: 'Euclidean', others: 'Maximum distance']
   
10. ###   lag        : time lag
    
11. ###
    
12. ###  References
    
13. ###  [1] Ma Jian, Sun Zengqi. Mutual information is copula entropy.
    
14. ###      arXiv:0808.0845, 2008.
    
15. ###  [2] Kraskov A, St\"ogbauer H, Grassberger P. Estimating mutual information.
    
16. ###      Physical review E, 2004, 69(6): 066138.
    
17. ###  [3] Ma Jian. Estimating Transfer Entropy via Copula Entropy.
    
18. ###      arXiv preprint arXiv:1910.04375, 2019.
    
19. ##################################################################################
    
20. copent<-function(x,k=3,dt=2){
    
21.   xc = construct_empirical_copula(x)
    
22.   -entknn(xc,k,dt)
    
23. }
    
24. 
    
25. construct_empirical_copula<-function(x){
    
26.   dimx = dim(as.matrix(x));
    
27.   rx = dimx[1];
    
28.   lx = dimx[2];
    
29.   
    
30.   xrank = x
    
31.   for(i in 1:lx){
    
32.     xrank[,i] = rank(x[,i])
    
33.   }
    
34.   
    
35.   xrank / rx
    
36. }
    
37. 
    
38. entknn<-function(x,k=3,dt=2){
    
39.   x = as.matrix(x)
    
40.   # get the dimension of data x
    
41.   N = dim(x)[1];
    
42.   d = dim(x)[2];
    
43.   
    
44.   g1 = digamma(N) - digamma(k);
    
45.   
    
46.   if (dt == 1){        # euciledean distance
    
47.     cd = pi^(d/2) / 2^d / gamma(1+d/2);       
    
48.     distx = as.matrix(dist(x));
    
49.   }
    
50.   else {        # maximum distance
    
51.     cd = 1;
    
52.     distx = as.matrix(dist(x,method = "maximum"));
    
53.   }
    
54.   
    
55.   logd = 0;
    
56.   for(i in 1:N){
    
57.     distx[i,] = sort(distx[i,]);
    
58.     logd = logd + log( 2 * distx[i,k+1] ) * d / N;
    
59.   }
    
60.   
    
61.   g1 + log(cd) + logd
    
62. }
    
63. 
    
64. ci<-function(x,y,z,k=3,dt=2){
    
65.   xyz = cbind(x,y,z)
    
66.   xz  = cbind(x,z)
    
67.   yz  = cbind(y,z)
    
68.   
    
69.   copent(xyz,k,dt) - copent(xz,k,dt) - copent(yz,k,dt)
    
70. }
    
71. 
    
72. transent<-function(x,y,lag=1,k=3,dt=2){
    
73.   xl = length(x)
    
74.   yl = length(y)
    
75.   if (xl > yl){
    
76.     l = yl
    
77.   }
    
78.   else {
    
79.     l = xl
    
80.   }
    
81.   x1 = x[1:(l-lag)]
    
82.   x2 = x[(lag+1):l]
    
83.   y1 = y[1:(l-lag)]
    
84. 
    
85.   ci(x2,y1,x1,k,dt)
    
86. }

复制代码