複数のマーカーについて一括してMakeWTable()を用いる

sphere MakeWTableS 分割表
  • 予定パッケージ名
  • 関数名
    • MakeWTableS
  • タイトル
MakeWTableS
  • 説明
複数のマーカーについて一括してMakeWTable()を用いる
  • 使用例
MakeWTableS(g,k)
  • ソース
MakeWTableS<-function(t,k=c(0,1)){
	n<-length(k)
	m<-matrix(0,length(t[1,]),length(t[,1])*n)
	for(i in 1:length(m[,1])){
		m[i,]<-MakeWTable(t[,i],k)
	}
	m
}
  • Rdファイル
\name{MakeWTableS}
\alias{MMakeWTableS}
%- Also NEED an '\alias' for EACH other topic documented here.
\title{
MakeWTable
}
\description{
複数のマーカーについて一括してMakeWTable()を用いる
}
\usage{
MakeWTableS(g,k)
}
\arguments{
\item{g}{サンプル数の行数、マーカー数の列数の行列}
\item{k}{カテゴリの重みづけベクトル}
}
%- maybe also 'usage' for other objects documented here.
\details{
%%  ~~ If necessary, more details than the description above ~~
}
\value{
マーカー数の行数、カテゴリ数ｘサンプル数の列数であるような行列
}
\references{
%% ~put references to the literature/web site here ~
}
\author{
%%  ~~who you are~~
}
\note{
%%  ~~further notes~~
}

%% ~Make other sections like Warning with \section{Warning }{....} ~

\seealso{
}
\examples{
Ns<-20
Nm<-5
k<-c(0,1)
num1<-sample(1:Ns,Nm,replace=TRUE)
num1<-sort(num1)
g<-matrix(0,Ns,Nm)
for(i in 1:Nm){
	g[,i]<-c(rep(0,Ns-num1[i]),rep(1,num1[i]))
}

MakeWTableS(t=g,k)

全サンプルを区別した分割表のための重みづけ表を作る

sphere MakeWTable 分割表
  • 予定パッケージ名
  • 関数名
    • MakeWTable
  • タイトル
MakeWTable
  • 説明
サンプルに値が与えられているときに、サンプル数ｘカテゴリ数の表に重みづけをつけることとする。
カテゴリには0:(カテゴリ数-1)の順序を付けるものとして重みづけする
  • 使用例
MakeWTable(g,k)
  • ソース
MakeWTable<-function(t,k=c(0,1)){
	n<-length(k)
	m<-matrix(0,n,length(t))
	for(i in 1:n){
		m[i,]<-t*k[i]
	}
	return(c(t(m)))
}
  • Rdファイル
\name{MakeWTable}
\alias{MMakeWTable}
%- Also NEED an '\alias' for EACH other topic documented here.
\title{
MakeWTable
}
\description{
サンプルに値が与えられているときに、サンプル数ｘカテゴリ数の表に重みづけをつけることとする。
カテゴリには0:(カテゴリ数-1)の順序を付けるものとして重みづけする
}
\usage{
MakeWTable(g,k)
}
\arguments{
\item{g}{サンプル数の長さのベクトル}
\item{k}{カテゴリの重みづけベクトル}
}
%- maybe also 'usage' for other objects documented here.
\details{
%%  ~~ If necessary, more details than the description above ~~
}
\value{
カテゴリ数ｘサンプル数であるようなベクトル
}
\references{
%% ~put references to the literature/web site here ~
}
\author{
%%  ~~who you are~~
}
\note{
%%  ~~further notes~~
}

%% ~Make other sections like Warning with \section{Warning }{....} ~

\seealso{
}
\examples{
Ns<-20
k<-c(0,1)
g<-sample(c(0,1,2),Ns,replace=TRUE)
w<-MakeWTable(g,k)

全サンプルを区別した分割表を作る

sphere MakeOTable 分割表
  • 予定パッケージ名
  • 関数名
    • MakeOTable
  • タイトル
MakeOTable
  • 説明
カテゴリ数がkでサンプル数がNのとき、個々のサンプルを１カラムと見立てた分割表は、表の要素が0か1かの表であり、列に関する和がすべて1であるような分割表である。
カテゴリ情報のベクトルとカテゴリのベクトルから、このような表を作る
  • 使用例
MakeOTable(p,k)
  • ソース
MakeOTable<-function(p,k=c(0,1)){
	n<-length(k)
	m<-matrix(0,n,length(p))
	for(i in 1:n){
		m[i,which(p==k[i])]<-1
	}
	m
}
  • Rdファイル
\name{MakeOTable}
\alias{MakeOTable}
%- Also NEED an '\alias' for EACH other topic documented here.
\title{
MakeOTable
}
\description{
カテゴリ数がkでサンプル数がNのとき、個々のサンプルを１カラムと見立てた分割表は、表の要素が0か1かの表であり、列に関する和がすべて1であるような分割表である。
カテゴリ情報のベクトルとカテゴリのベクトルから、このような表を作る
}
\usage{
MakeOTable(p,k)
}
\arguments{
\item{p}{サンプル数の長さのベクトル。要素はカテゴリ値}
\item{k}{カテゴリのベクトル}
}
%- maybe also 'usage' for other objects documented here.
\details{
%%  ~~ If necessary, more details than the description above ~~
}
\value{
行数がカテゴリ数、列数がサンプル数であるような行列
}
\references{
%% ~put references to the literature/web site here ~
}
\author{
%%  ~~who you are~~
}
\note{
%%  ~~further notes~~
}

%% ~Make other sections like Warning with \section{Warning }{....} ~

\seealso{
}
\examples{
Ns<-20
k<-c(0,1,3,"a")
p<-sample(k,Ns,replace=TRUE)
MakeOTable(p,k)

モンテカルロで自由度１の多重検定検出力

sphere MoteCarloPowerChisq 多次元球 分割表
  • Sweaveを用いて作成した解説メモはこちら
  • 予定パッケージ名
  • 関数名
    • SpherePower
  • タイトル
SpherePower
  • 説明
真値座標を指定、その指定倍距離の点を真値とする。
検定閾値距離を指定。
自由度１テストの法線ベクトルのセットを指定。
その最大値を統計量とする、多重検定を実施しているシチュエーション

そのうえで、モンテカルロ法にて、検定閾値を超える統計量を観察する。
  • 使用例
SPout<-SpherePower(df=df,Kv=Kv,Ks=Ks,CheckStats=CheckStats,Tiis=Tiis,nperm=nperm)
  • ソース
SpherePower<-function(df,Kv,Ks,CheckDists,Tiis,nperm=1000){
	Ks[which(Ks==0)]<-0.0001 # Ks=0は0から少しずらしておく。計算の都合による。
	# 同一方向、指定長さのベクトルのセットを作る
	Kss<-matrix(rep(Kv,length(Ks)),nrow=length(Ks),byrow=TRUE)*Ks
	# テストベクトルの数
	Ntest<-length(Tiis[,1])
	# 乱方向ベクトルのセットを発生する
	rss<-RandomSphere(df=df,n=nperm)
	# 出力格納ベクトル
	# Ksの値数ｘ検定統計量の値数のパワー値を算出する
	PowerOut<-matrix(0,length(Ks),length(CheckDists))
	# Ksの値ごとに計算を繰り返す
	for(ik in 1:length(Ks)){
		Pii<-Kss[ik,]
		# nperm個の方向について、指定検定統計量より大なる統計量を得る確率を格納する
		AccumProbs <- matrix(0, nperm, length(CheckDists))
		# 乱方向ごとに計算する
		for (i in 1:nperm) {
			rs <- matrix(rss[i,],nrow=1)
			tmptheta <- acos(cosxy(t(Pii), rs))
			C <- nrm(Pii) * cos(tmptheta)
			D <- nrm(Pii) * sin(tmptheta)
			As <- t(Tiis %*% t(rs))
			Bs <- t(Tiis %*% Pii)
			current <- which(Bs == max(Bs))[1]
			SerialTopT <- c(current)
			currentR <- 0
			SerialR <- c(currentR)
			SerialStat <- c(Bs[current])
			CheckList <- 1:Ntest
			CheckList <- CheckList[which(CheckList != current)]
			while (length(CheckList) > 0) {
				tmpks <- -(Bs[current] - Bs[CheckList])/(As[current] - 
				As[CheckList])
				tmpList <- which(tmpks > currentR)
				if (length(tmpList) > 0) {
					current <- CheckList[which(tmpks == min(tmpks[tmpList]))][1] # 次の等高線を作る線分のID
					currentR <- min(tmpks[tmpList]) # 対立仮説点からの距離
					SerialTopT <- c(SerialTopT, current) # IDをリストに
					SerialR <- c(SerialR, currentR) # 対立仮説点からの距離をリストに
					SerialStat <- c(SerialStat, As[current] * currentR + Bs[current]) # 統計量のリスト
				}	
				CheckList <- CheckList[tmpList] # 検討するべきIDのリストを更新
				CheckList <- CheckList[which(CheckList != current)] # 今登録したIDを除く
			}
			# intMyX:交叉すれば1
			# intMyX3:外から内へ入る場合に-1,内から外へ出る場合に+1
			# 線分は半開区間で、その始点は含み、終点は含まない
			intMyX<-intMyX3 <- matrix(0,length(CheckDists),length(SerialTopT))
			AnswerPs<-rep(0,length(CheckDists))
			for(j in 1:length(CheckDists)){
				tmpSerialStat<-SerialStat-CheckDists[j]
				#tmpSerialStat[which(abs(tmpSerialStat)<10^(-10))]<-0
				#print(tmpSerialStat)
				if(tmpSerialStat[1]>=0){
					AnswerPs[j]<-1
				}
#print("pre")
#print(AnswerPs)
				if(length(tmpSerialStat)>1){
					for(jj in 1:(length(tmpSerialStat)-1)){
						if(tmpSerialStat[jj]==0){
							intMyX[j,jj]<-1
							intMyX3[j,jj]<-sign(As[SerialTopT[jj]])
						}else{
							if(tmpSerialStat[jj]*tmpSerialStat[jj+1]<0){
								intMyX[j,jj]<-1
								intMyX3[j,jj]<-sign(As[SerialTopT[jj]])
							}
						}
					}
				}
				
				if(tmpSerialStat[length(tmpSerialStat)]==0){
					intMyX[j,length(tmpSerialStat)]<-1
					intMyX3[j,length(tmpSerialStat)]<-sign(As[SerialTopT[length(tmpSerialStat)]])
				}else{
					if(tmpSerialStat[length(tmpSerialStat)]*As[SerialTopT[length(tmpSerialStat)]]<0){
						intMyX[j,length(tmpSerialStat)]<-1
						intMyX3[j,length(tmpSerialStat)]<-sign(As[SerialTopT[length(tmpSerialStat)]])
					}
				}
			}


			zeros <- which(apply(abs(intMyX), 1, sum) == 0) # 交叉しないCheckDistsID？
			Crossings <- t(outer(CheckDists, Bs[SerialTopT], FUN = "-"))/As[SerialTopT]
			#print(intMyX3)
			#print(intMyX)
			#print(zeros)
			ChiProbs <- pchisq(Crossings^2, df, lower.tail = FALSE)
			#print(Crossings^2)
			#print(ChiProbs)
			AnswerPs <- AnswerPs + apply(t(intMyX3) * ChiProbs, 2, sum)
			AnswerPs[zeros]<-AnswerPs[zeros] <- (sign(As[SerialTopT[length(SerialTopT)]]) + 1)/2
			#print(AnswerPs)
			AnswerPs[which(AnswerPs<0)]<-1+AnswerPs[which(AnswerPs<0)]
			AnswerPs[which(AnswerPs>1)]<-AnswerPs[which(AnswerPs>1)]-1
			AccumProbs[i, ] <- AnswerPs
		}
		PowOut <- apply(AccumProbs, 2, sum)/nperm
		PowerOut[ik,]<-PowOut
	}

	list(df=df,Kv=Kv,Ks=Ks,CheckDists=CheckDists,Tiis=Tiis,nperm=nperm,p.out=PowerOut)
}
  • Rdファイル
\name{SpherePower}
\alias{SpherePower}
%- Also NEED an '\alias' for EACH other topic documented here.
\title{
SpherePower
}
\description{
真値座標を指定、その指定倍距離の点を真値とする。
検定閾値距離を指定。
自由度１テストの法線ベクトルのセットを指定。
その最大値を統計量とする、多重検定を実施しているシチュエーション

そのうえで、モンテカルロ法にて、検定閾値を超える統計量を観察する。
}
\usage{
SPout<-SpherePower(df=df,Kv=Kv,Ks=Ks,CheckStats=CheckStats,Tiis=Tiis,nperm=nperm)
}
\arguments{
\item{df}{自由度・次元}
\item{Kv}{真値方向ベクトル、単位ベクトルなら解釈が容易}
\item{Ks}{真値の原点からの距離}
\item{CheckDists}{検定の統計量(原点からの距離)。この二乗がカイ二乗値}
\item{Tiist}{法線ベクトルを行ベクトルとする行列}
\item{nperm}{モンテカルロ乱点数}
}
%- maybe also 'usage' for other objects documented here.
\details{
%%  ~~ If necessary, more details than the description above ~~
}
\value{
\item{df}{自由度・次元}
\item{Kv}{真値方向ベクトル、単位ベクトルなら解釈が容易}
\item{Ks}{真値の原点からの距離}
\item{CheckDists}{検定の統計量(原点からの距離)。この二乗がカイ二乗値}
\item{Tiist}{法線ベクトルを行ベクトルとする行列}
\item{nperm}{モンテカルロ乱点数}
\item{p.out}{Ks値を行に、CheckDistsを列として、CheckDists半径球の外側に点を得る確率}
}
\references{
%% ~put references to the literature/web site here ~
}
\author{
%%  ~~who you are~~
}
\note{
%%  ~~further notes~~
}

%% ~Make other sections like Warning with \section{Warning }{....} ~

\seealso{
}
\examples{
df<-sample(1:10,1)
t<-runif(1)*2*pi
#Kv<-c(cos(t),sin(t),0)
Kv<-runif(df)
Kv<-Kv/sqrt(sum(Kv^2))
Ks<-sample(1:50,1)
CheckDists<-seq(from=0,to=sqrt(Ks^2+df-1)*1.2,length.out=50)
#CheckDists<-Ks+df-1

#Tiis<-matrix(c(1,0,0,0,1,0,0,0,1),3,3,byrow=TRUE)
Nt<-1
Nt<-sample(1:1000,1)
library(MCMCpack)
Tiis<-rdirichlet(Nt,rep(1,df))
Tiis<-rbind(Tiis,-Tiis)
nperm<-100
spout<-SpherePower(df=df,Kv=Kv,Ks=Ks,CheckDists=CheckDists,Tiis=Tiis,nperm=nperm)

plot(c(spout$p.out),type="l")
par(new=TRUE)
plot(pchisq(CheckDists^2,df,Ks^2,lower.tail=FALSE),col=2,type="l")

[カイ二乗検定モンテカルロでカイ二乗検定の検出力

sphere MoteCarloPowerChisq 多次元球 分割表
  • 予定パッケージ名
  • 関数名
    • MoteCarloPowerChisq
  • タイトル
MoteCarloPowerChisq
  • 説明
指定次元、真値としてカイ二乗統計量を与え、棄却水準(p値もしくは、カイ二乗値)に対応するパワーを算出
非心カイ二乗分布での算出と、モンテカルロ算出との２方法を返す
  • 使用例
MoteCarloPowerChisq(K=15,df=4)
  • ソース
MoteCarloPowerChisq<-function(as=c(0.1,0.01,0.001),K,df,N=10000,p=TRUE){
	if(p){
		as<-qchisq(as,df,lower.tail=FALSE)
	}
	mr<-K
	xs<-RandomSphere(df,n=N) # 球面上の均一ランダム点座標を作る
	C<-xs[,1] # そのコサイン
	S<-sqrt(1-C^2) # そのサイン
	# 非心カイ二乗分布で出すパワーと球面乱数で出すパワーを格納
	p1<-p2<-rep(0,length(a))

for(xx in 1:length(as)){
	a<-as[xx]
	tmp<-a^2-mr^2*S^2
	tmp[which(tmp<0)]<-0
	tmp2<-mr*C
	tmp3<-sqrt(tmp)
	tmp4<-tmp2-tmp3
	tmp5<-tmp2+tmp3
	tmp4<-abs(tmp4)
	tmp5<-abs(tmp5)
	R12<-matrix(c(tmp4,tmp5),ncol=2)
	R1<-apply(R12,1,min)
	R2<-apply(R12,1,max)
	#plot(R1,R2)
# 検定閾値と真値との大小で処理が異なる
	if(a<mr){
		Q2<-0.5 + (pchisq(R1^2,df)+pchisq(R2^2,df,lower.tail=FALSE))/2
	}else{
		Q2<-(pchisq(R1^2,df,lower.tail=FALSE)+pchisq(R2^2,df,lower.tail=FALSE))/2
	}

p1[xx]<-mean(Q2)
p2[xx]<-pchisq(a^2,df,mr^2,lower.tail=FALSE)

}
# 一致する
list(p.ncchisq=p1,p.sim=p2)
}
  • Rdファイル
\name{MoteCarloPowerChisq}
\alias{MoteCarloPowerChisq}
%- Also NEED an '\alias' for EACH other topic documented here.
\title{
MoteCarloPowerChisq
}
\description{
指定次元、真値としてカイ二乗統計量を与え、棄却水準(p値もしくは、カイ二乗値)に対応するパワーを算出
非心カイ二乗分布での算出と、モンテカルロ算出との２方法を返す
}
\usage{
MoteCarloPowerChisq(K=15,df=4)
}
\arguments{
\item{as}{検定有意水準のp値リスト、もしくは、カイ二乗値リスト}
\item{K}{カイ二乗値の真値}
\item{df}{自由度}
\item{N}{モンテカルロ乱数点の個数}
\item{p}{TRUEはasがp値、FALSEはasがカイ二乗値}
}
%- maybe also 'usage' for other objects documented here.
\details{
%%  ~~ If necessary, more details than the description above ~~
}
\value{
\item{p.ncchisq}{非心カイ二乗分布によるパワー}
\item{p.sim}{モンテカルロシミュレーションによるパワー}
}
\references{
%% ~put references to the literature/web site here ~
}
\author{
%%  ~~who you are~~
}
\note{
%%  ~~further notes~~
}

%% ~Make other sections like Warning with \section{Warning }{....} ~

\seealso{
}
\examples{
MoteCarloPowerChisq(K=3,df=4)
#####
Ntrial<-100
Ma<-30
as<-seq(from=0,to=Ma,length.out=Ntrial) # 検定の有意水準はa^2できることにする
N<-10000 # 乱数点の数指定
k<-sample(2:10,1) # 次元を指定
mr<-runif(1)*5 # その原点からの距離
xs<-RandomSphere(k,n=N) # 球面上の均一ランダム点座標を作る
C<-xs[,1] # そのコサイン
S<-sqrt(1-C^2) # そのサイン
# 非心カイ二乗分布で出すパワーと球面乱数で出すパワーを格納
p1<-p2<-rep(0,Ntrial)

for(xx in 1:Ntrial){
	a<-as[xx]
	tmp<-a^2-mr^2*S^2
	tmp[which(tmp<0)]<-0
	tmp2<-mr*C
	tmp3<-sqrt(tmp)
	tmp4<-tmp2-tmp3
	tmp5<-tmp2+tmp3
	tmp4<-abs(tmp4)
	tmp5<-abs(tmp5)
	R12<-matrix(c(tmp4,tmp5),ncol=2)
	R1<-apply(R12,1,min)
	R2<-apply(R12,1,max)
	plot(R1,R2)
# 検定閾値と真値との大小で処理が異なる
	if(a<mr){
		Q2<-0.5 + (pchisq(R1^2,k)+pchisq(R2^2,k,lower.tail=FALSE))/2
	}else{
		Q2<-(pchisq(R1^2,k,lower.tail=FALSE)+pchisq(R2^2,k,lower.tail=FALSE))/2
	}

p1[xx]<-mean(Q2)
p2[xx]<-pchisq(a^2,k,mr^2,lower.tail=FALSE)

}
# 一致する
plot(p1,p2)

球面上の均一分布の乱数発生

sphere RandomSphere 多次元球 分割表 カイ二乗検定
  • 予定パッケージ名
  • 関数名
    • RandomSphere
  • タイトル
RandomSphere
  • 説明
任意次元の、指定半径多次元球面上の均一分布の発生
  • 使用例
RandomSphere(6)
  • ソース
RandomSphere<-
function (df = 3, r = 1, n = 100) 
{
    rs <- matrix(rnorm(df * n), nrow = n)
    rs/sqrt(apply(rs^2, 1, sum)) * r
}
  • Rdファイル
\name{RandomSphere}
\alias{RandomSphere}
%- Also NEED an '\alias' for EACH other topic documented here.
\title{
RandomSphere
}
\description{
任意次元の、指定半径多次元球面上の均一分布の発生
}
\usage{
RandomSphere(4)
}
\arguments{
df = 3, r = 1, n = 100
\item{df}{自然数　次元}
\item{r}{半径}
\item{n}{乱数点の数}
}
%- maybe also 'usage' for other objects documented here.
\details{
%%  ~~ If necessary, more details than the description above ~~
}
\value{
点の数の行ｘ次元の列の行列
}
\references{
%% ~put references to the literature/web site here ~
}
\author{
%%  ~~who you are~~
}
\note{
%%  ~~further notes~~
}

%% ~Make other sections like Warning with \section{Warning }{....} ~

\seealso{
}
\examples{
df<-3
x<-RandomSphere(df)
plot(as.data.frame(x))
}

自由度１の重みづけベクトルを自由度次元ベクトルに変換する

sphere TestDf1 分割表 カイ二乗検定
  • 予定パッケージ名
  • 関数名
    • TestDf1
  • タイトル
TestDf1
  • 説明
分割表の自由度１検定とは、分割表のセルに重みづけ係数をつけることである。
その長さNxMの重みづけベクトルを、標準化正単体座標系(自由度次元)の法線ベクトルに変換する
  • 使用例
dfCoordinate(O,Obs=TRUE)
  • ソース
TestDf1<-
function (O, W, Sp = NULL, simple = TRUE) 
{
    if (is.null(Sp)) {
        Sp <- RegularSpherize(O)
    }
    preT <- W
    score <- sum(preT * Sp$D)
    T <- apply(c(t(preT)) * Sp$X, 2, sum)
    Ti <- Sp$Ui %*% T
    Tii <- Ti/sqrt(Sp$EigenOut$values)
    Tiist <- Tii/sqrt(sum(Tii^2))
    cost <- 0
    if (sum(Sp$Pii^2) != 0) {
        cost <- sum(c(Tiist) * c(Sp$Pii))/sqrt(sum(Sp$Pii^2))
    }
    if (simple) {
        return(list(Tiist = Tiist, cost = cost, Qii = Tiist * 
            cost * sqrt(sum(Sp$Pii^2))))
    }
    else {
        K1 <- Sp$K * cost^2
        p <- pchisq(K1, 1, lower.tail = FALSE)
        list(Sp = Sp, Score = score, Tiist = Tiist, K1 = K1, 
            cost = cost, statistic = K1, p.value = p, Qii = Tiist * 
                cost * sqrt(sum(Sp$Pii^2)))
    }
}
  • Rdファイル
\name{TestDf1}
\alias{TestDf1}
%- Also NEED an '\alias' for EACH other topic documented here.
\title{
TestDf1
}
\description{
分割表の自由度１検定とは、分割表のセルに重みづけ係数をつけることである。
その長さNxMの重みづけベクトルを、標準化正単体座標系(自由度次元)の法線ベクトルに変換する
}
\usage{
TestDf1(O,W)
}
\arguments{
\item{O}{NxM 観察テーブル行列}
\item{W}{NxM 重みづけ行列}
\item{Sp}{RegularSpherize(O)が予め計算してあれば、それを与えてもよい}
\item{simple}{TRUE/FALSE,FALSEであれば、法線ベクトルのみを返し、TRUEであれば、検定結果も返す}
}
%- maybe also 'usage' for other objects documented here.
\details{
%%  ~~ If necessary, more details than the description above ~~
}
\value{
\item{Sp}{RegularSpherize(O)の出力}
\item{Score}{統計量(正負どちらもありえる)}
\item{Tiist}{標準化座標でのテスト法線ベクトル}
\item{K1}{標準化座標を用いて算出した統計量(非負)}
\item{cost}{標準化座標での観察テーブルベクトルとテスト法線ベクトルとのなす余弦}
\item{statistic}{統計量(非負)}
\item{p.value}{自由度１検定p値}
\item{Qii}{観察表ベクトルをテスト法線ベクトルへの射影ベクトル}
}
\references{
%% ~put references to the literature/web site here ~
}
\author{
%%  ~~who you are~~
}
\note{
%%  ~~further notes~~
}

%% ~Make other sections like Warning with \section{Warning }{....} ~

\seealso{
}
\examples{
O<-matrix(runif(6),2,3)
W<-matrix(runif(6),2,3)
TestDf1(O,W,simple=FALSE)
}