R에서 BayesLiNGAM

R에서 BayesLiNGAM



07May'20: Written

소개



현재 인과 검색 접근법에 사용되는 LiNGAM은 R과 Python으로 구현됩니다.

다양한 LiNGAM 모델의 R and/or Python
htps : // / s. 오, ぇ. 코 m / 해 / 시시 미즈 06 / ㄹ m

이번에는 그 중에서도 BayesLiNGAM을 이용한 모델을 R로 사용해 본다.
원리는 다음 문장에 기술되어 있다.

htps : // pdfs. 세만치 cs 쵸 r. rg / 7122 / 03d0696 fb5f03b1d046 6db59f7 4087cf. pdf
htps //w w. js 가짜. jst. . jp/아리치 cぇ/pj사이/J사이 2014/0/J사이 2014_2G13/_pdf

도입 환경



OS: Mojave (version; 10.14.6)
R 버전: 3.6.3 (2020-02-29)
R Studio 버전: 1.2.5033

절차


  • 패키지 다운로드 및 압축 해제
  • 소스 코드 재 작성
  • 실행에 필요한 다른 패키지 설치
  • 동작 확인

  • 패키지 다운로드 및 압축 해제



    아래 사이트의 "Code Package"에서 "Version 1.1"을 다운로드 (07May'20).
    htps //w w. cs. 헤르신키. 후 / g 로우 p / 네우 로인 f / ㄱ ㄴ m / 바 ぇ s ㎃ m /

    "bayeslingam.tar.gz"가 다운로드되기 때문에 이것을 압축 해제.
    → "bayeslingam"

    소스 코드 재작성



    이대로 R로 기동하려고 해도 Error가 되므로, 일부 소스 코드를 재기록한다.
    다시 쓰는 것은 2개.
    · loud.R
    · loadbayeslingam.R

    loud.R



    bayeslingam/main/loud.R의 일부를 재작성.
    원래는 이하.
    cauzality_path<<-"./../.."
    

    수정 후는 이하.
    cauzality_path<<-"/Users/(作業ディレクトリ)/bayeslingam"
    

    작업 디렉토리에는 선택 사항입니다.

    loadbayeslingam.R



    bayeslingam/bayeslingam/loadbayeslingam.R의 일부를 다시 씁니다.
    원래는 이하.
    common_dir<-sprintf("%s/trunk/common",cauzality_path)
    
    bayeslingam_dir<-sprintf("%s/trunk/bayeslingam",cauzality_path)    
    
    rbayeslingam_dir<-sprintf("%s/trunk/bayeslingam/rbayeslingam",cauzality_path)
    
    rlingam_dir<-sprintf("%s/trunk/rlingam",cauzality_path)
    
    rdeal_dir<-sprintf("%s/trunk/rdeal",cauzality_path)
    
    rpc_dir<-sprintf("%s/trunk/rpc",cauzality_path)
    
    etudag_dir<-sprintf("%s/trunk/ETUDAG",cauzality_path)
    

    trunk/를 지웁니다.
    수정 후는 이하.
    common_dir<-sprintf("%s/common",cauzality_path)
    
    bayeslingam_dir<-sprintf("%s/bayeslingam",cauzality_path)    
    
    rbayeslingam_dir<-sprintf("%s/bayeslingam/rbayeslingam",cauzality_path)
    
    rlingam_dir<-sprintf("%s/rlingam",cauzality_path)
    
    rdeal_dir<-sprintf("%s/rdeal",cauzality_path)
    
    rpc_dir<-sprintf("%s/rpc",cauzality_path)
    
    etudag_dir<-sprintf("%s/ETUDAG",cauzality_path)
    

    실행에 필요한 다른 패키지 설치



    다음 패키지를 설치합니다.
    이번에는 'fastICA'만으로도 가능.


    Package
    Feature


    deal
    GH-algorithm

    mclust
    MoG Bayeslingam

    combinat
    PC-algorithm

    ggm
    PC-algorithm

    fastICA
    Standard LiNGAM

    gtools
    Some method's in ETUDAG use this package (not used by BL)

    install.packages('PACKAGE')
    library(PACKAGE)
    

    동작 확인



    RStudio에서 bayeslingam을로드합니다.
    setwd(~/作業ディレクトリ)
    library('fastICA')
    source('bayeslingam/main/loud.R')
    loud()
    

    시도해보십시오.
    이번에는 이하의 인과 그래프를 추정해 본다.



    우선은 데이터 생성.
    x1 <- runif(5000, 0, 1)
    x2 <- x1 + runif(5000, 0, 1)
    x3 <- x1 + x2 + runif(5000, 0, 1)
    x4 <- x1 + x3 + runif(5000, 0, 1)
    

    데이터 프레임에 저장하고 BayesLiNGAM을 실시.
    d <- data.frame(x1=x1, x2=x2, x3=x3, x4=x4)
    result <- bayeslingam(d)
    

    DAGs를 살펴보십시오.
    > result$DAGs
           [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
      [1,]    1    2    3    4    0    0    0    0    0     0
      [2,]    1    2    3    4    0    0    0    0    0     1
      [3,]    1    2    3    4    0    0    0    0    1     0
      [4,]    1    2    3    4    0    0    0    0    1     1
      [5,]    1    2    3    4    0    0    0    1    0     0
    

    1열부터 4열째: 각 변수의 정렬 순서.
    e.g.) [2,]는 x1, x2, x3, x4를 나타낸다.

    5열부터 10열째: 유향선의 유무를 0,1로 나타낸다. 1열 내지 4열로 배열된 변수로부터 왼쪽 변수로의 유향선의 유무를 나타낸다.
    e.g.) [5,]에서 [,5]는 x1에서 x2로의 유향선 없음. [,6]는 x1에서 x3까지의 유향선 없음. [,7]은 x1에서 x4까지의 유향선 없음. [,8]는 x2에서 x3까지의 유향선 있다. [,9]는 x2에서 x4까지의 유향선 없음. [,10]은 x3에서 x4까지의 유향선 없음.

    또한, result$prob는 그 행이 나타내는 인과 그래프인 추정 확률을 나타낸다.

    마지막으로, 인과 그래프를 만든다.
    library(igraph)
    par(mfrow=c(2,3))
    prob <- round(result$prob,digits=4) * 100
    n.node <- max(result$components$node)
    case <- nrow(result$DAGs)
    node <- result$DAGs[,1:n.node]
    res <- as.matrix(result$DAGs[,-c(1:n.node)],nrow=case)
    name <-paste("X",1:n.node,sep="")
    
    for(i in order(prob,decreasing=T)[1:6]){
        amat <- matrix(0,n.node,n.node)
        index <- node[i,]
        amat[lower.tri(amat,diag=FALSE)] <- res[i,]
        amat <- t(amat[index,index])
        g <- graph.adjacency(amat,mode="directed")
        E(g)$label <- NA
        pos <- layout.circle(g)
        rot <- matrix(c(0,1,-1,0),nrow=2,byrow=T)
        la <- pos %*% rot
        if(n.node == 2)la <- matrix(c(-1,0.5,1,0.5),,nrow=2,byrow=T)
        plot(g,layout=la,vertex.size=70,vertex.label=name)
        mtext(paste(prob[i],"%"),line=0)
    }
    

    이 결과는 다음과 같습니다.



    98%에서 인과 그래프를 확인할 수 있었다.

    결론



    기본적으로 @m__k의 "LiNGAM의 베이지안 접근법 인 BayesLiNGAM을 움직여 보았다"(https://qiita.com/m__k/items/288910248efafa106863)을 따라 실행하면 좋았습니다.
    그리고, 인자가 5개 이상의 경우에 「greedybayeslingam()」를 사용하는 것 같지만 잘 되지 않았다.

    참조



    ・LiNGAM의 베이즈적 어프로치인 BayesLiNGAM을 움직여 보았다 ( @m__k 씨)
    htps : // 이 m/m__k/있어 ms/288910248
    · NagoyaR19_인과 분석 LiNGAM을 시험해 보았다 (요시노 아츠시 씨)
    h tps://d ゔぇ. 오, ぇ. 코 m / 후 / d / 0B5 및 qfxcbRc WZ 또는 T Ql 호아이 1F 1 jRWJ 후 k9wSmRTdVF · ぃ에 w
    ・A→B인가, B→A인가를 데이터로부터 간파할 수 있을까? (LiNGAM 시뮬레이션을 해 보았다)
    https://blog.albert2005.co.jp/2015/02/24/a%E2%86%92b%E3%81%AA%E3%81%AE%E3%81%8B%E3%80%81b% E2%86%92a%E3%81%AA%E3%81%AE%E3%81%8B%E3%82%92%E3%83%87%E3%83%BC%E3%82%BF%E3% 81%8B%E3%82%89%E8%A6%8B%E6%8A%9C%E3%81%8F%E3%81%93%E3%81%A8%E3%81%AF/
    · BayesLiNGAM (IMOTO Yusuke 씨)
    htps : // / s. 오, ぇ. 이 m/ゔぃえw/yー이모토/%에 5% 82% 99%에 5% BF% 98%에 9% 8C% B2/바에 s ぃ가 m

    좋은 웹페이지 즐겨찾기