★ R -- U 検定（base で用意されてはいるが，練習で） ★

　355　R -- U 検定（base で用意されてはいるが，練習で）　　青木繁伸　　2002/01/24 (木) 23:37
　　357　Re: R -- U 検定（base で用意されてはいるが，練習で）　　青木繁伸　　2002/01/25 (金) 00:00
　　　359　Re^2: R -- U 検定（base で用意されてはいるが，練習で）　　sb812109　　2002/01/25 (金) 00:06
　　　　363　Re^3: R -- U 検定（base で用意されてはいるが，練習で）　　青木繁伸　　2002/01/25 (金) 10:46

355.　R -- U 検定（base で用意されてはいるが，練習で）　　青木繁伸　　2002/01/24 (木) 23:37 
u.test <- function(x, y)
{
    n <- (n1 <- length(x))+(n2 <- length(y))
    u <- n1*n2+n1*(n1+1)/2-sum(rank(z <- c(x, y))[1:n1])
    u <- min(u1, n1*n2-u1)
    ties <- table(z)
    z <- abs(u-(e <- n1*n2/2))/sqrt(v <- n1*n2*(n^3-n-sum(ties^3-ties))/(12*(n*n-n)))
    ans <- c(u, e, v, z, pnorm(z, low=F)*2)
    names(ans) <- c("U", "E[U]", "Var[U]", "Z", "P")
    ans
}

# 例
u.test(c(1.2, 1.5, 1.8, 2.6), c(1.3, 1.9, 2.9, 3.1, 3.9))
#         U       E[U]     Var[U]          Z          P 
# 4.0000000 10.0000000 16.6666667  1.4696938  0.1416447

　　　　　[このページのトップへ]

357.　Re: R -- U 検定（base で用意されてはいるが，練習で）　　青木繁伸　　2002/01/25 (金) 00:00 
インターフェースのバリエーション

# data データ行列
# d    測定値のある列番号
# f    群識別値のある列番号
# g1   第一群を表す値
# g2   第二群を表す値

u.test2 <- function(data, d, f, g1, g2)
{
    x <- data[,d][data[,f] == g1]
    y <- data[,d][data[,f] == g2]
    u.test(x,y)
}

# 使用例（1列目:測定値，2列目:群番号）
# ファイルから読んでも良い

data <- matrix(c(
1.2 1
1.5 1
1.8 1
2.6 1
1.3 2
1.9 2
2.9 2
3.1 2
3.9 2
), ncol=2, byrow=TRUE)

u.test2(data, 1, 2, 1, 2)

　　　　　[このページのトップへ]

359.　Re^2: R -- U 検定（base で用意されてはいるが，練習で）　　sb812109　　2002/01/25 (金) 00:06

沢山の R スクリプト有難うございます。Library(Aoki)として利用
させて戴いています。

ところで，246 の kappa 統計量ですが，次の様なメッセージが返っ
てきます。
Error in pnorm(z, lower = F) : Object "z" not found

宜しくお願いします。
(管理者注:もとの発言を書庫へ移すときに修正してあります）

　　　　　[このページのトップへ]

363.　Re^3: R -- U 検定（base で用意されてはいるが，練習で）　　青木繁伸　　2002/01/25 (金) 10:46
> 沢山の R スクリプト有難うございます。Library(Aoki)として利用
> させて戴いています。

恐縮です。

> ところで，246 の kappa 統計量ですが，次の様なメッセージが返っ
> てきます。
> Error in pnorm(z, lower = F) : Object "z" not found


短くしようとして，ミスしてしまいました。以下のように修正させていただきます。

(管理者注:書庫へ入れるときに，もとの発言も修正してあります） 

    sk0 <- sqrt(pe/(n*qe))
    result <- c(kappa, sk, sk0, kappa-qnorm(0.975)*sk, kappa+qnorm(0.975)*sk, kappa/sk0, pnorm(z,lower=F))

を

    sk0 <- sqrt(pe/(n*qe))
    z<-kappa/sk0
    result <- c(kappa, sk, sk0, kappa-qnorm(0.975)*sk, kappa+qnorm(0.975)*sk, z, pnorm(z,lower=F))
に

　　　　　[このページのトップへ]

● 「統計学関連なんでもあり」の過去ログ--- 017 の目次へジャンプ
● 「統計学関連なんでもあり」の目次へジャンプ
● 直前のページへ戻る