将后验分布叠加在卡明的猫眼可视化均值上

Superimpose posterior distribution on mean like cat's eye visualization from cumming

我个人喜欢卡明的猫眼可视化,它将采样分布叠加在点估计上:

我还想用 Scripts of Kruschke (2015) 获得的后验分布来做这个:

plotMCMC( mcmcCoda , data=myData , 
      compValMu=100.0 , ropeMu=c(99.0,101.0) ,
      compValSigma=15.0 , ropeSigma=c(14,16) ,
      compValEff=0.0 , ropeEff=c(-0.1,0.1) ,
      pairsPlot=TRUE , showCurve=FALSE ,
      saveName=fileNameRoot , saveType=graphFileType )

 # Set up window and layout:


openGraph(width=6.0,height=8.0*3/5)
  layout( matrix( c(2,3,5, 1,4,6) , nrow=3, byrow=FALSE ) )
  par( mar=c(3.5,3.5,2.5,0.5) , mgp=c(2.25,0.7,0) )
  # Select thinned steps in chain for plotting of posterior predictive     curves:
  nCurvesToPlot = 20
  stepIdxVec = seq( 1 , chainLength , floor(chainLength/nCurvesToPlot) )
  # Compute limits for plots of data with posterior pred. distributions
  y = data
  xLim = c( min(y)-0.1*(max(y)-min(y)) , max(y)+0.1*(max(y)-min(y)) )
  xBreaks = seq( xLim[1] , xLim[2] , 
                 length=ceiling((xLim[2]-xLim[1])/(sd(y)/4)) )
  histInfo = hist(y,breaks=xBreaks,plot=FALSE)
  yMax = 1.2 * max( histInfo$density )
  xVec = seq( xLim[1] , xLim[2] , length=501 )
  #-----------------------------------------------------------------------------
  # Plot data y and smattering of posterior predictive curves:
  histInfo = hist( y , prob=TRUE , xlim=xLim , ylim=c(0,yMax) , breaks=xBreaks,
                   col="red2" , border="white" , xlab="y" , ylab="" , 
                   yaxt="n" , cex.lab=1.5 , main="Data w. Post. Pred." )
  for ( stepIdx in 1:length(stepIdxVec) ) {
    lines(xVec, dt( (xVec-mu[stepIdxVec[stepIdx]])/sigma[stepIdxVec[stepIdx]], 
                    df=nu[stepIdxVec[stepIdx]] )/sigma[stepIdxVec[stepIdx]] , 
          type="l" , col="skyblue" , lwd=1 )
  }
  text( max(xVec) , yMax , bquote(N==.(length(y))) , adj=c(1.1,1.1) )
  #-----------------------------------------------------------------------------
  histInfo = plotPost( mu , cex.lab = 1.75 , showCurve=showCurve ,
                       compVal=compValMu , ROPE=ropeMu ,
                       xlab=bquote(mu) , main=paste("Mean") , 
                       col="skyblue" )
  #-----------------------------------------------------------------------------
  histInfo = plotPost( sigma , cex.lab=1.75 , showCurve=showCurve ,
                       compVal=compValSigma , ROPE=ropeSigma , cenTend="mode" ,
                       xlab=bquote(sigma) , main=paste("Scale") , 
                       col="skyblue" )
  #-----------------------------------------------------------------------------
  effectSize = ( mu - compValMu ) / sigma
  histInfo = plotPost( effectSize , compVal=compValEff ,  ROPE=ropeEff ,
                       showCurve=showCurve , cenTend="mode" ,
                       xlab=bquote( ( mu - .(compValMu) ) / sigma ),
                       cex.lab=1.75 , main="Effect Size" ,
                       col="skyblue" )
  #-----------------------------------------------------------------------------  
  postInfo = plotPost( log10(nu) , col="skyblue" , # breaks=30 ,
                       showCurve=showCurve ,
                       xlab=bquote("log10("*nu*")") , cex.lab = 1.75 , 
                       cenTend="mode" ,
                       main="Normality" ) #  (<0.7 suggests kurtosis)

最后看起来像这样的平均值

是否可以像卡明那样将直方图叠加在均值上?

我已经为另一个应用程序创建了一个名为 plotViolin 的函数,它应该可以执行您想要的操作,或者至少可以根据您的应用程序的具体情况轻松修改。您可以在此处链接的 R 脚本的一半位置找到该函数:https://osf.io/wt4vf/