abikoushi

library(MCMCpack)
library(exact2x2)
library(ggplot2)
library(patchwork)

X <-matrix(c(12,5,6,12), nrow=2)
fpfun <- function(phi){exact2x2::exact2x2(x=X, or=exp(phi))$p.value}
rs <- rowSums(X)
cs <- colSums(X)
flfun <- function(phi){MCMCpack::dnoncenhypergeom(x = X[1,1], cs[1],cs[2],rs[1], exp(phi))}

abikoushi

library(ggplot2)
set.seed(123)
x <- rweibull(50,1.5) - log(2)^(1/1.5)
qweibull(0.5,1.5)- log(2)^(1/1.5) #中央値が0の分布
sx <- jitter(sort(x)) #ちょうど0をなくすため
pv_w <- sapply(sx,function(m)wilcox.test(x, mu=m)$p.value)
pv_s <- sapply(sx,function(m)binom.test(sum(x>m), length(x))$p.value)
df4p <- data.frame(x=x,sx=sx,sign=pv_s,wilcox=pv_w)

m <- median(x)

abikoushi

X <- as.matrix(iris[,-5])
d <- dist(X)
h1 <- hclust(d^2, method = "ward.D")
h2 <- hclust(d, method = "ward.D2")
h1$height <- sqrt(h1$height)

#png("dendro2.png")
par(mfrow=c(1,2))
plot(h1)
plot(h2)

abikoushi

Y<-iris[,-5]
svdout <- svd(log(Y))
cname <- hcl.colors(3, alpha = 0.1)
#png("iris_svd.png")
matplot(t(svdout$u),
        type = "l",lty=1,
        col = cname[iris$Species],
        xlab = "PC", ylab = "value")
legend("topleft", legend = levels(iris$Species),
       col=hcl.colors(3), lty=1)

abikoushi

#include <RcppArmadillo.h>
using namespace Rcpp;
using namespace arma;
// [[Rcpp::depends(RcppArmadillo)]]


int rcate0(const int & K){
  int x = 0;
  double U = R::runif(0,1);
  if(U > 1/K){

abikoushi

#https://okumuralab.org/~okumura/stat/2by2.html

library(MCMCpack)

confint_fisher <- function(X, level=0.95){
  rs <- rowSums(X)
  cs <- colSums(X)
  alpha <- 1-level
  testfun_up <- function(phi){
    pall <- MCMCpack::dnoncenhypergeom(x = NA, cs[1],cs[2],rs[1], exp(phi))

abikoushi

library(ggplot2)
library(dplyr)

doVB <- function(y,X,a=0.5,b=0.01,iter=100){
  K <- ncol(X)
  sxy = t(X)%*%y
  lambda <- rgamma(K,1,1)
  ahat <- rep(a,K)
  bhat <- rep(b,K)
  for (i in 1:iter) {

abikoushi

library(ggplot2)
library(dplyr)


wilcox_p <- function(mu,x,prob=0.5){
  x2 <- x-mu
  cmb <- combn(10,2)
  u <-apply(cmb,2,function(i)(x2[i[1]]+x2[i[2]])/2)
  u <- c(u,x2)
  pos <- sum(u>0)

abikoushi

#include <RcppArmadillo.h>
using namespace Rcpp;
using namespace arma;
// [[Rcpp::depends(RcppArmadillo)]]

// [[Rcpp::export]]
mat rcate(const int N, const rowvec & p){
  int K = p.n_cols;
  rowvec cump = cumsum(p);
  mat x(N, K);

abikoushi

doVB <- function(y,X,a=0.5,b=0.001,iter=200){
  K <- ncol(X)
  sxy = t(X)%*%y
  lambda <- rgamma(K,1,1)
  ahat <- rep(a,K)
  bhat <- rep(b,K)
  for (i in 1:iter) {
    for(j in 1:K){
      ahat[j] = sxy[j]+a
      bhat[j] = X[,j]%*%exp(X[,-j,drop=FALSE]%*%log(lambda[-j]))+b