スピアマンの順位相関係数の関数化

質問

タイトル通り、スピアマンの順位相関係数の関数化ができずに困っています。

統計ソフトのＲを使っています。本来ならcor(a,b)で0.02857143とでます。これを関数化せよ、ということです。

自分で以下のようにやってみましたが上手くいきません。どこが違って、どのようにしたらいいか、アドバイスをお願いします。

> a <- c(3,2,1,5,4,6)
> b <- c(2,6,3,5,1,4)
> n <- length(a)

> test <- function(a,b,n){ 1-(6*sum(a-b)^2)/
(n*(n+1)*(n-1))}

回答No.2

grothendieck

> test
とだけ入力するとtestの定義を表示することになります。testを実行するためには関数に渡すオブジェクトを指定する必要があります。下に実際にRで実行した結果を示します。

> x <- c(3,2,1,5,4,6)
> y <- c(2,6,3,5,1,4)
> k <- length(x)
> test <- function(a,b,n) {1-6*sum((a-b)^2)/(n*(n+1)*(n-1))}
> test(x,y,k)
[1] 0.02857143

function(a,b,n)でのa,b,nは関数を定義するための変数なので、関数を実行する時にa,b,nに代入されるオブジェクトを指定します。

お礼

オブジェクトを指定しないとだめなんですね！勉強になりました。>test(a,b,n)としたところ、きちんと実行されました。ありがとうございました！

回答No.1

grothendieck

sum(a-b)^2 とすると sum(a-b) の方が先に計算されて
　sum(a-b)^2 = (sum(a-b))^2 = 0
になるようです。(a-b)^2 を先に計算するために
　sum((a-b)^2)
としたらうまく行きました。なお、下記URLにはRのスピアマン相関係数の計算には欠測値の処理に問題があること、およびその対策があります。

参考URL：

http://aoki2.si.gunma-u.ac.jp/R/spearman.html

補足

ありがとうございます＾＾

上の書き方で上手く動いたんですか？私は上のように書いて、> test とやっても結果が表示されずにtest<-function・・・　が出てしまいます。表示の仕方が悪いのでしょうか・・？