【R言語】統計解析フリーソフトＲ 第6章【GNU R】 [無断転載禁止]©2ch.net
■ このスレッドは過去ログ倉庫に格納されています
>>309
CRANはここ
https://cran.r-project.org
プロジェクトの方はESETの誤検知っぽいんだよな >>310
してねぇよｗｗ
totalvirusで調べると一件引っかかるし、Redirector検知だから、しょうがないね col=rgb(runif(1),runif(1),runif(1),runif(1))で色指定すると
走らせるたびに色がちがっておもしろい。
hist(rnorm(100),col=rgb(runif(1),runif(1),runif(1),runif(1))) stanやJAGSのコードでgamma関数を使おうとして
y = gamma(x)
と、書いたらエラーになった。
stanだと y=tgamma(x)、JAGSだとy=exp(loggam(x))で動作した。 ある大学の入学者男女の比率は１であるという帰無仮説を検定する課題が花子と太郎に課された。
花子は５０人を調査できたら終了として入学者を５０人をみつけて１８人が女子であるという結果を得た。
帰無仮説のもとで
５０人中１８人が女子である確率は 0.01603475
これ以下になるのは５０人中０〜１８人と３２〜５０人が女子の場合なので
両側検定して
> sum(dbinom(c(0:18,32:50),50,0.5))
[1] 0.06490865
> binom.test(18,50,0.5)$p.value
[1] 0.06490865
で帰無仮説は棄却できないと結論した。
http://i.imgur.com/XDIp9rM.png
一方、本番と十八番が好きな太郎は一人ずつ調べて１８人めの女子がみつかったところで調査を終えることにした。
１８人めがみつかったのは花子と同じく５０人めであった。
帰無仮説のもとで
１８人がみつかるのが５０人めである確率は0.005772512
これ以下になるのは２３人以下５０人以上番めで女子１８人めがみつかった場合なので
両側検定して
pnb=dnbinom(0:999,18,0.5)
> 1 - sum(pnb[-which(pnb<=dnbinom(50-18,18,0.5))]) # < 0.05
[1] 0.02750309
http://i.imgur.com/K3T7utr.png
で帰無仮説は棄却される。
どちらの検定が正しいか、どちらも正しくないか？
検定する意図によってｐ値が変わるのは頻度主義統計の欠陥といえるか？
花子の横軸は裏口入学者数、太郎の横軸はサンプル数なので
サンプルでの裏口入学率を横軸にして９５％信頼区間を示す。
花子の検定での信頼区間は0.36〜0.72で18/50を含む、p=0.06491
http://i.imgur.com/SeTLk8K.jpg
太郎の検定での信頼区間は0.375〜0.72で18/50を含まない、p= 0.0275
http://i.imgur.com/tNzlfxe.jpg
主観である、検定の中止の基準の差でｐ値や信頼区間が変化するのは変だという批判である。 （訂正）
ある大学の入学者男女の比率は１であるという帰無仮説を検定する課題が花子と太郎に課された。
花子は５０人を調査できたら終了として入学者を５０人をみつけて１８人が女子であるという結果を得た。
帰無仮説のもとで
５０人中１８人が女子である確率は 0.01603475
これ以下になるのは５０人中０〜１８人と３２〜５０人が女子の場合なので
両側検定して
> sum(dbinom(c(0:18,32:50),50,0.5))
[1] 0.06490865
> binom.test(18,50,0.5)$p.value
[1] 0.06490865
で帰無仮説は棄却できないと結論した。
http://i.imgur.com/XDIp9rM.png
一方、十八という数字が好きな太郎は一人ずつ調べて１８人めの女子がみつかったところで調査を終えることにした。
１８人めがみつかったのは花子と同じく５０人めであった。
帰無仮説のもとで
１８人がみつかるのが５０人めである確率は0.005772512
これ以下になるのは２３人以下５０人以上番めで女子１８人めがみつかった場合なので
両側検定して
pnb=dnbinom(0:999,18,0.5)
> 1 - sum(pnb[-which(pnb<=dnbinom(50-18,18,0.5))]) # < 0.05
[1] 0.02750309
http://i.imgur.com/K3T7utr.png
で帰無仮説は棄却される。
どちらの検定が正しいか、どちらも正しくないか？
検定する意図によってｐ値が変わるのは頻度主義統計の欠陥といえるか？
花子の横軸は女子数、太郎の横軸はサンプル数なので
サンプルでの女子の割合を横軸にして９５％信頼区間を示す。
花子の検定での信頼区間は0.36〜0.72で18/50を含む、p=0.06491
http://i.imgur.com/SeTLk8K.jpg
太郎の検定での信頼区間は0.375〜0.72で18/50を含まない、p= 0.0275
http://i.imgur.com/tNzlfxe.jpg
主観である、検定の中止の基準の差でｐ値や信頼区間が変化するのは変だという批判である。 kainokousiki<-function(a,b,c){return (-b+sqrt(b^2-4*a*c))/(2*a)} #解の公式
kainokousiki(1,-5,6)
でrunすると3じゃなくて6を返すんだけど、どこが間違ってる？ 機械学習をきっかけにPythonに逆転された感じだね 『Rを使った〜』だとPythonじゃないのかよって思うよね ずっとRだけでPython触ったこと無いけど、覚え直す価値ある？
環境構築からもう面倒なイメージ Python自体は:と直後のインデントさえ気を付ければかなり簡単
3系は数が全て小数扱いなので楽
Anacondaというパッケージでインストールすれば、今流行りのJupyter Notebookという開発環境で対話的にコーディングできる（Rも使える） アナコンダてのがRぽくできるのね、ありがとう
dplyrやggplot2みたいに素人でも簡単便利だといいんだけど
Pythonはオブジェクト志向ぽいしすぐ諦めそう >>335
誤爆？
何をブラウズする環境が無くなったの？ RStudio Severあるで。Dockerで使うのがいいんじゃないかな？！ 統計学とウェブ解析を交えて実践的な勉強と練習を
したいのですが、おすすめな書籍やサイトはありますか。
実際に解析ツールや分析ツールを用いて
自分で分析解析してから
解答を見て適切な手順や方法、考察を
解説してくれるものが良いです。
統計学は統計検定2級の知識はありますが
ウェブ解析はテキスト読んだだけです。 R studioは日本語コメント書く度にIMEが無効になったりカーソルがずれたり黒文字の予測変換が黒背景と重なって見えなくなったりと散々だわ 会社規定なんでしかたない
UNIX環境使えるのが羨ましい Docker for WindowsでRStudioサーバー動かせば？ >>342
IME無効になるのは俺だけじゃなかったと知ってほっとした。 ストアアプリも同じ症状でるからRStudio固有の問題でなくWindows環境の不治の病だと思ってる Run any R code you like. There are over three thousand R packages preloaded.
https://rdrr.io/snippets/ Fisher test検定時に
p<2.2e-16
と表示されるんですが、これより小さい値の指数桁数を正確に表記する方法教えて下さい。
例えば5.8e-35となるようにです。 fisher.test関数の返り値はlist型で、その中にp.valueという名前でp値が格納されているから$演算子を使って直接参照するか、broom::tidy関数に返り値を渡してdata.frame形式で出力してやれば見れる >>352
>>354が言う直接的な参照
> fisher.test(matrix(c(1,120,130,2),2))$p.value
[1] 1.691912e-69 352です。よく分かりました！
ありがとうございます！ >>356
技術的な助言をしたけど、学術的に言えば、
p < 0.01 は全て p < 0.01 として、具体的なp値を考える意味はないと思うよ。
一部の例外的な研究分野を除いて(e.g., 遺伝統計学)。 はい、まさにその例外的な分野で使おうとしてます。ありがとうございます。 Rのガンマ関数はいくつでオーバーフローするかやってみた。
> i=1
> while(gamma(i)!=Inf){
+ i=i+1
+ }
Warning message:
In gamma(i) : value out of range in 'gammafn'
> i
[1] 172
> gamma(172)
[1] Inf
Warning message:
value out of range in 'gammafn'
> gamma(171)
[1] 7.257416e+306 matplot()で折れ線グラフ描いたときに、X軸をカテゴリで示したいのですが、
可能でしょうか？
例えばtemp <- c("0時間","8時間","24時間","48時間")として、
matplot()の引数にtempをとるやり方です。
他にもやり方あれば教えてください。 >>371
matplot(..., xaxt="n")
axis(1, at=seq(along=temp), lab=temp) >>372
遅くなりましたがありがとうございました。
できました！ 特定の長方形の中に複数の長方形を最小面積で敷き詰める平面充填に関するパッケージってありませんかね # jonckheereテストを書いてみた
jonckheere <- function(L,
alternative = c("two.sided", "increasing", "decreasing"),
cat=TRUE){
# L : list of vectors A1,A2,...,Ak, with assumed tendency
How.Many.Greater.Pairs <- function(A,B){ # How many pairs of A[i] > B[j], count as 0.5 when equal,
n.a = length(A)
n.b = length(B)
how.many.greater.pairs = 0
for(i in 1:n.a){
for(j in 1:n.b){
how.many.greater.pairs = how.many.greater.pairs+ifelse(A[i]==B[j],0.5,A[i]>B[j])
}
}
return(how.many.greater.pairs)
}
Sum.of.Greater.Pairs <- function(L){ #L=list(A1,,,,Ak),A1 < A2 < A3,..,< Ak : vector
k = length(L)
comb = combn(1:k,2) # possible combinaition of pairs to compare
n.comb = ncol(comb) # how many combinations
J = 0 # sum of greater pairs
for(i in 1:n.comb){
J = J + How.Many.Greater.Pairs(L[[comb[1,i]]],L[[comb[2,i]]])
}
return(J)
}
J = Sum.of.Greater.Pairs(L)
n = sapply(L,length)
N = sum(n)
EJ = (N^2-sum(n^2))/4
VJ = (N^2*(2*N+3)-sum(n^2*(2*n+3)))/72
Z = (J-EJ)/sqrt(VJ)
alternative = match.arg(alternative)
p.value = switch(alternative, 'two.sided' = 2 * min(pnorm(Z), pnorm(-Z), 0.5),
'increasing' = pnorm(Z),
'decreasing' = pnorm(-Z))
if(cat){
cat( 'p.value = ', p.value,'\n')
cat('alternative hypothesis: ' ,alternative,'\n')
}
invisible(p.value)
} >>388
実態調査か何か？
<-と=は挙動が違う場合があるので、使い分けていますが、
代入はどっちかと問われたら、無論 <- または ->
なお、
> 1 -> x
これはエラーにならないけど、
> 1 = x
1 = x でエラー: 代入の左辺が不正 (do_set) です
これはエラー >>389
俺は基本=派。
関数の定義は
z.test <- function(x,n=16,sigma=1){
z=sqrt(n)*mean(x)/sigma
2*pnorm(abs(z),lower=FALSE)
}
と書いている。 ＞1 = x でエラー
当たり前
そんな使い方なんてするかよ
他言語と同じく=一文字の方がすっきりしてイイ こういうのが紛らわしいから、俺は = 推奨。
x <- 1
if(x <- 1) print('YES')
if(x < -1) print('YES') やったことなかったので関数の初期値設定に<-を使うとどうなるかやってみた。
まず、= の場合
> z.test <- function(x,n=16,sigma=1){
+ z=sqrt(n)*mean(x)/sigma
+ 2*pnorm(abs(z),lower=FALSE)
+ }
> z.test(1:3)
[1] 1.244192e-15
<- で初期値設定すると、エラー
> z0.test <- function(x,n<-16,sigma<-1){
Error: unexpected assignment in "z0.test <- function(x,n<-"
> z=sqrt(n)*mean(x)/sigma
Error in mean(x) : object 'x' not found
> 2*pnorm(abs(z),lower=FALSE)
Error in pnorm(abs(z), lower = FALSE) : object 'z' not found
> } 俺は見栄えがいいと思って関数定義には<-を使っているけど
= でも通常に動作する。
> z.test = function(x,n=16,sigma=1){
+ z=sqrt(n)*mean(x)/sigma
+ 2*pnorm(abs(z),lower=FALSE)
+ }
> z.test(1:3)
[1] 1.244192e-15
<- 推奨の人に聞きたいのだけど
<- でないと動作しないってことあるのだろうか？ RStudioつかってりゃ［Alt+-］で簡単入力
それと引数の指定は代入じゃないと思うんだが感覚が違うよかな？ >>396
essだとアンダースコアを入れるを勝手に「<-」になる。
>>394
それは代入ではなく、関数の規定値の設定。
規定値の設定は「=」と決められているので、「<-」はアウト。
ただし、関数を実行するときには「<-」を使うことができる。
> mean(x<-1:10)
[1] 5.5
> x
[1] 1 2 3 4 5 6 7 8 9 10 =に拘り続けて、<-に変える瞬間の君たちがたまらないよ >>398
ありがとうございます。
そんな使い方ができたのですね。
こんなのしか知りませんでした。
> f <- function() {
+ x<<-1:10
+ mean (x)
+ }
> f()
[1] 5.5
> x
[1] 1 2 3 4 5 6 7 8 9 10 一画面しかグラフがないのに
Hit <Return> to see next plot:
と出るの鬱陶しいな、と思っていたら
par(ask=FALSE)
と設定しておけばいいんだな。 同名のdataがあるのでパッケージを::で指定してもうまくいかなかった。
data(netmeta::parkinson)
でなくて
data(parkinson, package = 'netmeta')
とするんだな。 hist(rnorm(100),col=rgb(runif(1),runif(1),runif(1),runif(1)))
だと、1色だけど
hist(rnorm(100),col=sample(colours(),(sample(1:10,1))))
なら1〜10色で表示される。 > hist(rnorm(100),col=rgb(runif(1),runif(1),runif(1),runif(1)))
何色に増やしてもいいがな。
何をやりたいのかな？
hist(rnorm(5000),col=apply(matrix(runif(80),4), 2, function(x){rgb(x[1],x[2],x[3],x[4])})) :::てどういう時に使うんだろ？
ソース読みたくて
library(BayesFactor)
ttestBF_indepSample
では表示されなかったが、
library(BayesFactor)
BayesFactor:::ttestBF_indepSample
だと出てきた。 ■ このスレッドは過去ログ倉庫に格納されています
