ryamadaの遺伝学・遺伝統計学メモ このページをアンテナに追加 RSSフィード

数学・コンピュータ関連の姉妹ブログ『ryamadaのコンピュータ・数学メモ』
京都大学大学院医学研究科ゲノム医学センター統計遺伝学分野のWiki
講義・スライド
医学生物学と数学とプログラミングの三重学習を狙う学習ツール
駆け足で読む○○シリーズ
ぱらぱらめくるシリーズ
カシオの計算機
オンライン整数列大辞典

2016-05-21 2のd乗組み合わせデータ構造のための情報幾何

[][][][][][]Dually Flat Manifolds

  • 昨日、階層構造を持つ確率分布のための情報幾何について少し書いた
  • 特に、¥{0,1¥}^dなd元分割表のためのe-平坦、m-平坦なパラメタの取り方について、丁寧に確認してみたい
  • Rでパラメタ変換をしてみるのが、手っ取り早そうなので、それでやってみる
  • Rmd
---
title: "Dual flat parameterization for {0,1}^d"
author: "ryamada"
date: "Friday, May 20, 2016"
output: html_document
---

# Brief introduction

{0,1}の値を取るd個の変数$X_i;i=1,2,...$があるときに、その値セットの取り方($\mathbf{x} = (0,0,...,0),(1,0,...),...,(1,1,...,1)$)は
$D=2^d$通りある。

今、この$D$通りの取り方の生起確率ベクトルを$freq=(f_1,...,f_D)$とすれば、
$\sum f_i = 1$の制約があり、自由度は$D-1$である。


このような階層構造を持つ組み合わせについては、情報幾何的な観点から有用な、パラメタの取り方が2通り知られている。

It is know that there are tow useful ways to parameterize the frequency with the hierarchical structure.

この文書の目的は、頻度ベクトルを2通りのパラメタに置き換えたり、パラメタの値を指定して、頻度ベクトルを作ったりする道具を提供すること。

This document is written to offer utility functions to transform frequency vectors to two types of parameters and their inverse.




詳しい理論は ["Information Geometry on Hierarchy of Probability
Distributions" (IEEE TRANSACTIONS ON INFORMATION THEORY, VOL. 47, NO. 5, JULY 2001)](https://www.google.co.jp/url?sa=t&rct=j&q=&esrc=s&source=web&cd=1&cad=rja&uact=8&ved=0ahUKEwja3aaXwOXMAhXHXaYKHVSXCzAQFgghMAA&url=http%3A%2F%2Fieeexplore.ieee.org%2Fiel5%2F18%2F20133%2F00930911.pdf%3Farnumber%3D930911&usg=AFQjCNFd9avfcHWB3jmeVpfvnuWzCKo5rQ&sig2=OBOG6MXZgpPJvHd2JksQIQ) を参照のこと。

You can learn its theory in ["Information Geometry on Hierarchy of Probability
Distributions" (IEEE TRANSACTIONS ON INFORMATION THEORY, VOL. 47, NO. 5, JULY 2001)](https://www.google.co.jp/url?sa=t&rct=j&q=&esrc=s&source=web&cd=1&cad=rja&uact=8&ved=0ahUKEwja3aaXwOXMAhXHXaYKHVSXCzAQFgghMAA&url=http%3A%2F%2Fieeexplore.ieee.org%2Fiel5%2F18%2F20133%2F00930911.pdf%3Farnumber%3D930911&usg=AFQjCNFd9avfcHWB3jmeVpfvnuWzCKo5rQ&sig2=OBOG6MXZgpPJvHd2JksQIQ)

## e-flat parameterization

$\mathbf{\theta}$は長さ$D-1$のパラメタベクトル。

$\mathbf{\theta}$ is a parameter vector with $D-1$ elements.

$\mathbf{\theta} = ((\theta_1,\theta_2,...,\theta_d,),(\theta_{1,2},...,\theta_{d-1,d}),...,(\theta_{1,2,...,d}))$

$\mathbf{x}=(x_1,...,x_d)$というパターンを取る確率
$$
\log{p(\mathbf{x})}=(x_1,...,x_d))= \sum \theta_i x_i + \sum \theta_{i,j} x_i x_j + ... + \theta_{1,2,...,d} x_1 x_2 ... x_d - \psi
$$

$psi$の値は$\mathbf{\theta}$が決まると(自由度がもうないので)一意に決まる。

The value of $\psi$ is not free because $\mathbf{\theta}$ has $D-1=df$ elements.

$\theta$ は1要素の寄与、要素ペアの寄与、要素トリオの寄与…に対応づく。


$\theta$ represents single effect, pairwise effect, triowise... and so on.

## m-flat parameterization

$\mathbf{\eta}$は、長さ$D-1$のパラメタベクトル。

$\mathbf{\eta}$ is a parameter vector with $D-1$ elements.

$\mathbf{\eta} = ((\eta_1,\eta_2,...,\eta_d,),(\eta_{1,2},...,\eta_{d-1,d}),...,(\eta_{1,2,...,d}))$

$\eta$の各要素は、変数(単独・組み合わせの)期待値を表す。(多元分割表で言えば)周辺頻度にも対応する。

$\eta$ corresponds to the expected value of single/combination of variables, or in other words, marginal frequency of multi-way table.

$$
\eta_{i} = Prob(x_i=1)\\
\eta_{i,j} = Prob(x_i=x_j=1)\\
...
\eta_{1,2,...,d} = Prob(x_1=x_2=...=x_d=1)
$$

# R codes in the ryamada22/Ronlyryamada github-package

# Utility functions and their usage example
```{r}
#library(devtools)
#install_github("ryamada22/Ronlyryamada")
library(Ronlyryamada)
example(my.haplotypes)
```

# Vignettes on functions

## Enumeration of $2^d$ patterns
d個の変数、$2^d$の場合の列挙。
1行目から順番に、値が1である変数の数が、0,1,...,dと増えるような並び方になっている。

$2^d$ patterns for d variables.
From the top row to the bottom, number of variables with value 1 increases like 0,1,...,d.


```{}
my.haplotypes <- function(d){
  n <- 2^d
	X <- matrix(0,d,n)
	cnt <- 1
	for(i in 1:d){
		tmp <- combn(d,i)
		for(j in 1:length(tmp[1,])){
			X[tmp[,j],cnt+j] <- 1
		}
		cnt <- cnt + length(tmp[1,])
	}
	return(X)
}
```
```{r}
d <- 3
my.haplotypes(d)
```

## Matrices to convert between frequency, theta and eta.

```{}
my.freq2eta.mat <- function(d){
  X <- my.haplotypes(d)
	n <- 2^d
	M <- matrix(0,n,n)

	for(i in 1:n){
		M[,i] <- apply(X >= X[,i],2,prod)
	}
	return(t(M[2:n,2:n]))
}
my.theta2freq.mat <- function(d){
  K <- my.freq2eta.mat(d)
	n <- 2^d
	M <- matrix(0,n,n)
	M[2:n,2:n] <- K
	M[,1] <- -1
	return(M)
}
my.freq2theta.mat <- function(d){
  M <- my.theta2freq.mat(d)
	return(solve(M))
}
# eta2freq.mat = inverse of freq2eta.mat
```
```{r}
my.freq2eta.mat(d)
my.theta2freq.mat(d)
my.freq2theta.mat(d)
```

## Transformations among three parameterizations

"freq", "$\theta$", "$eta$"3通りのパラメタ表現。

"$\theta$"の場合、自動的に決まるD番目の値を計算しておく方がよい。ほかの2つのパラメタシステムの場合は足して1になるところからの差、と言う形で簡単に計算できるが、$psi$はそうではなく面倒くさいから。


"freq" stands for regular frequency vector with $2^d=D$ elements.

"$theta$" stands for e-flat parameters.

"$eta$" stands for m-flat parameters.

For $\theta$, the residual parameter "\psi" should be also calculated because it is a bit difficult than other two systems' D-th value, even though it is automatically determined. 

```{}
my.theta2psi <- function(theta){
  len <- length(theta)
	d <- log(len+1,2)
	M <- my.theta2freq.mat(d)
	#tmp <- M %*% c(0,theta)
	tmp <- my.Inf.matrix.mult(M,c(0,theta))
	delta <- sum(exp(tmp))
	psi <- log(delta)
	return(list(psi=psi,M=M))
}
#' @export

my.theta2freq <- function(theta){
	out <- my.theta2psi(theta)
	psi.theta <- c(out$psi,theta)
	M <- out$M
	return(exp(my.Inf.matrix.mult(M,psi.theta)))
}
#' @export

my.freq2thetapsi <- function(freq){
	d <- log(length(freq),2)
	M <- my.freq2theta.mat(d)
	#return(M %*% log(freq))
	return(my.Inf.matrix.mult(M,log(freq)))
}
#' @export

my.freq2eta <- function(freq){
	d <- log(length(freq),2)
	K <- my.freq2eta.mat(d)
	return(K %*% freq[-1])
}
#' @export

my.eta2freq <- function(eta){
	d <- log(length(eta)+1,2)
	K.inv <- solve(my.freq2eta.mat(d))
	tmp <- K.inv %*% eta
	return(c(1-sum(tmp),tmp))
}
#' @export

my.theta2eta <- function(theta){
	freq <- my.theta2freq(theta)
	return(my.freq2eta(freq))
}
#' @export

my.eta2thetapsi <- function(eta){
	freq <- my.eta2freq(eta)
	return(my.freq2thetapsi(freq))
}
```

## Infinite values の扱い

freqベクトルに0が含まれると、$\theta$が無限大になる。

無限大を含めて計算できるようにするために、少し、工夫が必要。そのユーティリティ関数。

無限大と非無限大とを分離して行列計算し、最後にまとめる。

```{r}
my.matrix.Inf.separator <- function(M){
  if(!is.matrix(M)){
		M <- matrix(M,ncol=1)
	}
	M.inf <- M.noninf <- matrix(0,length(M[,1]),length(M[1,]))
	infs <- which(is.infinite(M))
	noninfs <- which(!is.infinite(M))
	M.noninf[noninfs] <- M[noninfs]
	M.inf[infs] <- sign(M[infs])
	return(list(M.noninf=M.noninf,M.inf=M.inf))
}

my.Inf.matrix.mult <- function(M1,M2){
	out1 <- my.matrix.Inf.separator(M1)
	out2 <- my.matrix.Inf.separator(M2)
	
	tmp.n.n <- out1$M.noninf %*% out2$M.noninf
	tmp.i.i <- out1$M.inf %*% out2$M.inf
	tmp.n.i <- out1$M.noninf %*% out2$M.inf
	tmp.i.n <- out1$M.inf %*% out2$M.noninf
	
	tmp.any.i <- sign(tmp.i.i + tmp.n.i + tmp.i.n)
	non.zero <- which(tmp.any.i!=0)
	tmp.any.i[non.zero] <- tmp.any.i[non.zero] * Inf
	tmp <- tmp.n.n + tmp.any.i
	return(tmp)
}
```

2016-05-20 Dually Flat Manifolds

[][][][][][]Dually Flat Manifolds

  • 資料
  • 確率分布を点としてもつ多様体を統計多様体という
  • 確率分布をn個のパラメタで表すことにすると、その統計多様体はn次元多様体
  • 多様体には多様体固有の特徴があって、それはリーマン計量だったりするわけで、局所座標系を入れる・入れないに関わらずその多様体の広がり具合というのは決まっている
  • 統計多様体では、座標系を入れる。パラメタで分布を表現することと同じ
  • パラメタを入れて、多様体全体をパラメタを使った座標表示にすると、格子が張り付いたことになる
  • この格子は、うまく張り付けると、ねじれたり捩れたりしないで素直に多様体上を流れていく
  • そんなよい感じの座標表示を"flat"と呼ぶ
    • 具体的には、3次元以上の多様体のときに、任意の3パラメタの間に延び広がり方・曲り方に偏りがないようなもののこと
    • 2通りのflatがあって、e-flat, m-flatと言うものがあり、この2つのflatは相互に密接な関係がある
    • またe-,m-flatsを表す式もよく似ている
    • e-flatの定義:確率分布が¥mathbf{¥theta}で表されている
      • E¥[¥frac{¥partial^2}{¥partial ¥theta_i ¥partial ¥theta_j}¥log{p(x,¥mathbf{¥theta})} ¥frac{¥partial}{¥partial ¥theta_k}¥log{p(x,¥mathbf{¥theta})}¥]=0
    • m-flatの定義:確率分布が¥mathbf{¥eta}で表されている
      • E¥[¥frac{1}{p(x,¥mathbf{¥eta})}¥frac{¥partial^2}{¥partial ¥eta_i ¥partial ¥eta_j}p(x,¥mathbf{¥theta}) ¥frac{¥partial}{¥partial ¥eta_k}¥log{p(x,¥mathbf{¥eta})}¥]=0
  • e-flatならm-flat、m-flatならe-flat
  • e-,f-両方のflatなパラメタ表現があるとよい(のだが、そうそう簡単に見つからないかもしれない)
  • ありがたいことに…
    • 多くの確率分布のパラメタ表現は、e-flatになるような表現で書きだせる。また、そのような表現を持つ確率分布を指数分布族と呼ぶ。そのような表現は
      • p(x,¥mathbf{¥theta})=e^{¥sum ¥theta_i k_i(x)- ¥psi (¥mathbf{¥theta})}
        • ただし、e-flatな分布関数が、すべて指数分布族なわけではない
    • また、m-flatになるような確率分布族も知られている。混合分布族
      • p(x,¥mathbf{¥eta})=¥sum ¥eta_i q_i(x) + (1-¥sum ¥eta_i) q_0(x)
        • ただし、m-flatな関数が、すべて混合分布族なわけでもない
  • e-flatなパラメタ表現がわかっても、m-flat表現がわからないとつまらないし、その逆も然り
    • ありがたいことに、ルジャンドル変換と呼ばれるルールで見つけられる。以下の変換ルールを見ればわかるように、¥mathbf{¥theta}¥mathbf{¥eta}とは、要素ごとに対応関係がある
      • ¥eta_i = ¥frac{¥partial}{¥partial ¥theta_i} ¥psi(¥mathbf{¥theta})= E¥[k_i(x)¥]
      • ¥theta_i = ¥frac{¥partial}{¥partial ¥eta_i} ¥phi(¥mathbf{¥eta})
        • ¥phi(¥mathbf{¥eta}) = E¥[¥log{p(x,¥mathbf{¥eta})}¥]
  • 対数尤度関数での関係
    • 確率分布を2通りのflatなパラメタ表現をすること、それらの相互関係について書いてきた
    • 確率分布関数を尤度関数とみなすとき、尤度関数の代わりに対数尤度関数を用いることも多い。確かに、上記の記載でもあちこちで対数化した関数が用いられていた
    • 以下では、対数尤度関数¥log{p(x,¥mathbf{¥theta})},¥log{p(x,¥mathbf{¥eta})}としてやり、そのパラメタ方向偏微分を成分とするベクトルを考えると、¥theta,¥etaとが「正規直交基底関係」のように見える(実際は、¥mathbf{¥theta}のそれぞれの成分同士は互いに直交ではなく、¥mathbf{¥eta}も同様なので、e-flat側とm-flat側との対応がきれいだ、ということになる
    • <e_i,e^*_j> = ¥delta_{ij}
  • Dual flatなパラメタ表記を用いると、次が言える
    • ¥psi(¥mathbf{¥theta}) + ¥phi(¥mathbf{¥eta}) + ¥mathbf{¥theta} ¥cdot ¥mathbf{¥eta} = 0
    • この式は、p(x,¥mathbf{¥theta}),p(x,¥mathbf{¥eta}という、「同じ関数」「多様体上の同じ点」の異なるパラメタ表記に関すること
    • ここで、「異なる点」についてこれを用いる
    • D¥[p:p’¥]=¥psi(¥mathbf{¥theta}) + ¥phi(¥mathbf{¥eta}’) + ¥mathbf{¥theta} ¥cdot ¥mathbf{¥eta}’というものが定義できて、これはp=p’のときに0でそうでないときには正の値をとる
    • これがKL divergenceとなっている
    • D¥[p:p’¥] = E_{¥mathbf{¥theta}}¥[¥log{¥frac{p(x,¥mathbf{¥theta})}{p(x,¥mathbf{¥theta}’)}}¥]
  • 関数間の遠近関係
    • 関数は統計多様体上の点なので、関数の違いを多様体上の「距離」として測りたい
    • しかしながら、KL divergenceの定義は、「自身」とのdivergenceを0となるような定義式を使って、「他者」のdivergenceを正の値にしており、2つの関数のどちらを基準にするかで値が変わる。そのことを受けて「距離」ではなく"divergence"と呼ぶ
    • e-,m-flatsという2通りのパラメタ系の良さは次の通り
      • e-系で表しているとき、e-flatなので、ごちゃごちゃ考えずとも、e-系のパラメタを使って比較的簡単にdivergenceを計量することができる
      • m-系の場合も同じ
      • 今、3つの関数p,p',p''があったときに、
      • p,p'について、e-系表現がわかっていれば、両者のdivergenceが簡単に測れる
      • さらに、p',p''について、m-系表現がわかっていれば、両者のdivergenceも簡単に測れる
      • p,p''の間のdivergenceは、片やe-系、片やm-系なので、簡単に測れないけれど、ここで、e-系とm-系とが「直交」していることが役に立ち、両者が単純に足し合わせればよいことが示せる
  • 関数間の遠近を統計量にする
    • KL divergence(対数尤度の差)の2倍の値はカイ二乗分布に従うことが知られており、尤度比検定に持ち込める
  • Flatなsubmanifolds
    • e-,m-flatsのよいところは、それぞれの座標系で、関数空間全体がベクトル空間のようになっていることである
    • 今、例えばe-系で表された関数空間全体のうち、e-系の線形和で表せる部分空間(submanifold)を取り出したとする
    • これはe-系表現があるからe-flatだし、e-flatならばm-flatではあるけれど、m-系で線形和に表せるかどうかというとそうはならない
      • e-系とm-系とは、相互に「直交」関係にあり、e-のすべての要素のそれぞれをm-系の全要素が分担して表現する、というような関係になっている(なっていがち)なために、e-の方の部分をとっても、m-の方は「個々の要素の寄与分」は小さくなっても、mの要素数を減らすわけにはいかないことなどを念頭に置くとわかりやすい
  • うまく辻褄の合うe-系submanifoldシリーズとm-系submanifoldシリーズ
    • e-系とm-系とのパラメタには、対応関係があるので、それについて注意をしつつ、パラメタを除いていくと(パラメタの値を固定していくと)、入れ子になったe系submanifoldsシリーズとm-系submanifoldsシリーズができる。それは、直交性を利用して活用することができる

2016-05-19 90分で学ぶ医学・ライフサイエンスの統計学手法のいまどき

[][][][][][][][][]手法・目的の分類オーバービュー〜90分で学ぶ医学・ライフサイエンスの統計学手法のいまどき

  • Statistical Bioinformatics: For Biomedical and Life Science Researchers(の目次)
  • Quality Control of High-Throughput Data
  • Statistical Tests, Statistical Significance, Error Controlling
  • Classification/Clustering
    • Unsupervised Learning
    • Supervised Learning
  • Multidimensional/Highdimensional Data and Data visualization
  • Statistical Models and Inference
  • Random value-based approaches
    • Resampling
    • MCMC
  • Graph and Network
  • Experimental Designs

[][][]Multiple-Comparison Issue

  • When you test multiple times, you should not believe nominal p-values of individual tests. 検定を複数行った場合には、個々の検定のp値をそのまま使って解釈できない
p <- runif(10^5)
hist(p)
plot(sort(p),pch=20,cex=0.1)
alpha <- 0.05
abline(v=length(p)*alpha,col="red")
abline(h=alpha,col="blue")
  • Multiple tests (100 tests) 100個の検定
n <- 100
m <- 1000
ps <- matrix(runif(n*m),n,m)
min.ps <- apply(ps,2,min)
hist(min.ps)
plot(sort(min.ps),pch=20,cex=0.1)
alpha <- 0.05
abline(v=length(min.ps)*alpha,col="red")
abline(h=alpha,col="blue")
  • False Disovery Rate
    • q-value not p-value
    • Some tests should be "significant"; many genes should be significantly diffrently expressed between cancer and normal cells.
    • The smallest p-value should be small enough, but the 2nd, 3rd,... smallest could be not so small and they should be judged by different threshold.
n.null <- 80
n.alt <- 20
p <- runif(n.null)
p <- c(p,pchisq(rchisq(n.alt,df=1,ncp=10),df=1,lower.tail=FALSE))
p <- sort(p)
plot(p)
q <- p.adjust(p,"BH")
k <- max(which(q<0.05))
plot(p,pch=20)
abline(0,p[k]/k,col="red")

f:id:ryamada22:20160518140425j:image

plot(q,pch=20)

f:id:ryamada22:20160518140426j:image

    • Further reading : Large-scale inference; When you have many observations all together, use their overall distribution/behavior for judgement of each.

[]High Dimensionality

f:id:ryamada22:20160518124000p:image

  • Sparse
    • 10^5cells with 5 markers
      • 0.1 ¥times 10grid; one cubicle only has 1 cell in average. It is too sparse to estimate density in a regular way.
  • No center, no common individuals

f:id:ryamada22:20160518125810j:image

N <- 10^4
x <- runif(N)-0.5
y <- runif(N)-0.5

s <- which(x^2+y^2<0.5^2)

plot(x,y,pch=20,cex=1)
points(x[s],y[s],pch=20,cex=1,col="red")
n <- 1:10
r <- 0.5
vol.hypersphere <- pi^(n/2)/gamma(n/2+1) * r^n
plot(n,vol.hypersphere,pch=20,cex=3)

f:id:ryamada22:20160518150548j:image

[]small n large p

  • 100 samples x 25000 genes
  • You can predict perfectly when you are allowed to use explanatory variables as many as sample size.
n <- 10
m <- 10
p <- sample(0:1,n,replace=TRUE)
g <- matrix(sample(0:1,n*m,replace=TRUE),n,m)
p
g
lm.out <- lm(p~g)
p.est <- predict(lm.out)

plot(p,p.est)
print(cbind(p,p.est))
print(cbind(p,round(p.est)))

print(round(lm.out$coefficients))
> p
 [1] 1 1 0 0 0 1 1 0 0 0
> g
      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
 [1,]    0    1    1    0    1    1    0    1    1     0
 [2,]    1    0    0    0    0    1    1    0    0     1
 [3,]    0    0    0    0    1    1    1    0    0     1
 [4,]    0    1    1    0    1    1    1    0    0     1
 [5,]    0    1    0    0    0    0    1    0    0     1
 [6,]    1    1    1    0    0    0    1    1    0     0
 [7,]    1    0    0    0    0    0    1    0    0     0
 [8,]    1    0    1    0    0    0    1    1    0     0
 [9,]    0    1    1    1    1    0    0    0    0     0
[10,]    1    0    0    0    1    1    1    0    0     0
> print(cbind(p,p.est))
   p         p.est
1  1  1.000000e+00
2  1  1.000000e+00
3  0 -3.333333e-01
4  0  8.881784e-16
5  0  3.333333e-01
6  1  6.666667e-01
7  1  6.666667e-01
8  0  3.333333e-01
9  0  1.110223e-16
10 0  3.333333e-01
> print(cbind(p,round(p.est)))
   p  
1  1 1
2  1 1
3  0 0
4  0 0
5  0 0
6  1 1
7  1 1
8  0 0
9  0 0
10 0 0
> 
> print(round(lm.out$coefficients))
(Intercept)          g1          g2          g3          g4          g5          g6          g7          g8 
          1           1           0           0          -1          -1           0          -1           0 
         g9         g10 
         NA          NA 
  • When you have p >> n explanatory variabls, you would find many variable sets, each of which perfectly predict dependent variable.

[]Noisy High-Throughput Biological Data

  • Noisy because
    • Biological phenomena are "noisy" ~ heterogeneity is the important feature of biology.
    • Experiments have many factors that add noise to data.
    • Highthroughput systems realize "highthoughput" by sacrificing preciseness somehow.
  • Noises are informative
    • Large-scale data have information on "shape" of noise and it can be used to quality-control/adjust the noisy data to purify "non-noisy signals".
  • 1000 samples, SNP chip for 50K SNP genotyping.
  • Success rate of individual samples.
N <- 1000
n1 <- 650
n2 <- 300

a <- rbeta(n1,50,2)
b <- rbeta(n2,100,20)
x <- rbeta(N-n1-n2,30,50)
hist(c(a,b,x),breaks = seq(from=0,to=1,length=40))

f:id:ryamada22:20160518145751j:image

[][]検定からベイズ推定へ Frequentist approach ~ tests and Bayesian approach

[][]MCMC

  • What is new?
  • This approach requires computers, that was not available in 20th century; that is why MCMC is being used for tasks that could not be solved by the methods before and also MCMC is being used for tasks that have been answered by the methods as well.
  • Link

[][][]時系列・空間とkrigingとか Time-series/Space and krigign etc.

  • Conventional analysis depends on independency among values
  • Time series/ spacial data have values that are next to each other are tightly associated and this feature should be absolutely considered.
  • At a glance

2016-05-18 Ricci flow と アルファ接続

[][][]Ricci flowとアルファ接続

  • ここ数日、統計多様体とその2つの平坦なパラメタの取り方について書いている
  • リーマン多様体として形を考えるときにRicci flowというのがあった
  • そこにもアルファ接続が出てくる
  • どういう関係なのか調べよう
  • こちらに資料

2016-05-17 情報幾何 Affine接続 捩れ双対接続 双対平坦

[][][][][]情報幾何 Affine接続 捩れ双対接続 双対平坦

  • 資料はこちら
  • リーマン多様体(多様体に計量が乗っている)がある
  • そこに滑らかにつながるベクトル場がある
  • ベクトル場を多様体上で微分したい
  • ベクトル場の微分をするとは、「あるベクトル場」を「別のベクトル場が定める方向」について微分してやり、「新たなベクトル場」を作ることだとする
  • したがって、ベクトル場の微分を考えるとき、同一の多様体上の3つのベクトル場が登場し、そこに首尾一貫する演算ルールを定めたい
  • XをY方向に微分するか、YをX方向に微分するか、で違いが出るか、出ないで収まるか、とか、そういう話になってくる
  • この微分のようなものは「ベクトル場xベクトル場→ベクトル場」という作業になる
  • そんな「ベクトル場xベクトル場→ベクトル場」の演算に線形性があるか、任意の多様体上関数に関して、よい扱いができるか、などを基準にしてにして、この演算の性質を定めたとき、それを「Affine 接続」と呼ぶ
  • Affine接続により、ベクトル場が多様体上を連なっていることになる。その連なり具合は、「曲がって」もよいけれど、「捩れてはいけない」という制約をいれる
  • そうすれば、2階の微分・変化はないことにできて、1階の関係だけにできるから(1階の関係だけ、と言っても、自由度・次元がnなら、そのトリオでn^3の関係を担保しないといけないので大変なのだが)
  • そんな曲りはするが捩れないAffine接続が、アルファ接続と呼ばれ、これは1変数 アルファを使って統一的に書けるもの
  • そのうち大事なのは、¥alpha=0¥alpha=¥pm 1の3つ
  • 0の場合は、Levi-Civita接続という名前を持ち、1のときはe-接続、-1のときはm-接続という名前を持つ
  • さらに、e-接続でありながら、その接続を表す、n^3の要素のすべての関係を表す係数が0のとき、それを「平坦」と呼び、e-平坦
  • 同様にm-接続でありながら、組み合わせ項が0になるものをm-平坦と呼ぶ(このとき、接続は捩れていないだけでなく、曲がってもいない)
  • また、このe-接続、m-接続は、n個の要素から3つを取り出す方法だが、微分幾何を思い出せば、3つを取り出すことと、その「裏側」にn-3個を取り出しつつ、ホッジ対応を指せることが可能なので、e-接続、m-接続はそれぞれ、別物であって、しかも大事
  • 別物であるときに、平坦性が自動的に担保されるかと言えばそうではないので、「平坦性」にも2通りあることになるが、それが両方とも平坦なものはよくて、そういうものが統計多様体で現れるので、それもよい
  • また、このe-接続、m-接続は、相互に補い合う関係(¥alpha,-¥alphaの関係)なのだが、それは、順序の入れ替えで符号が変わることの対応のようなもの
  • さらに¥alpha=0の場合は、その符号変換もない「対称性」を備えた、「一番きれいな」接続で、だからLevi-Civita接続は、一意
  • また、順序が変わることの影響を受ける¥pm 1のe-,m-接続は、それぞれ、「順序が変わることから、『距離』の定義を満たさないが、直交は満足している(ピタゴラスの定理を使えるような)divergenceがあり、それが実はKL-divergence」