110842
統計学関連なんでもあり

統計学に関する質問とか意見などなんでも書き込んでください。回答があるか,回答が正しいかは保証できません。
レポートや宿題の答えを求める記事や,その他不適切と判断した記事は,管理者の独断で即刻削除します。
ローマ数字,丸付き数字などのフォントセット依存文字および半角カタカナは使わないでください。
記事を引用する際には,適切に引用してください(全文引用はしないでください)。
問題が解決した(結局解決しなかった)場合は,その旨報告していただければ,コメントした人や他の読者に参考になるでしょう。


[トップに戻る] [利用上の注意] [ワード検索] [過去ログ] [統計学自習ノート] [管理用]
latest article: No. 22328, 2017/04/23(Sun) 12:46
おなまえ
タイトル
コメント
URL
添付ファイル
暗証キー (英数字で8文字以内)

自分の投稿記事を修正・削除する場合 --> 処理 記事No 暗証キー


二値データと連続変数の相関関係について
  投稿者:成瀬 2017/04/21(Fri) 16:11 No. 22325
初心者です。
基本的な質問で申し訳ございません。
二値データと連続変数の相関関係をみるには、二値データをダミー変数に置き換えて数値化して解析すると教わったのですが、これは正しい方法なのでしょうか?
二値データを数値化したところで、数値に意味がないわけなので、正しい方法とは思えないと考えています。

二値データと連続変数の相関関係をみる方法をご教示お願い申し上げます。

Re: 二値データと連続変数の相関関係について
  投稿者:課長 2017/04/22(Sat) 10:14 No. 22326
仰るとおり、1,0に意味は無いかと思います。計算上便利なだけでしょう。しかし、0,1を取るデータの数は意味があるのではないでしょうか?実際の計算では0や1を取る値の数、そして連続データの平均などの値を使うと思います。点双列相関係数やスピアマンの相関係数でいいんじゃないでしょうか?ちょっと自身がありませんので他の方の回答を期待します。

Re: 二値データと連続変数の相関関係について
  投稿者:青木繁伸 2017/04/22(Sat) 21:14 No. 22327
正しい方法です。
名前もあります。
点双列相関係数といいます。

しかし,それは,簡便化した計算方法を用いるものとして,「特別に?」定義づけられたものに過ぎず,実体は計算方法も含め「ピアソンの積率相関係数」にほかなりません。

「名義尺度を数値データとみて計算する」のは,うわべだけ見れば,「そんな馬鹿なことないよなあ」と思うかも知れませんが,名義尺度も二値データの場合は特別です。0/1 に限らず 1/2 であろうと50/1000 であろうと,標準化してしまえば同じ(標準化というのが,文字通り標準化ですから)。

相関係数を計算するには,間隔尺度,比尺度であること,というのがよく知られていて,名義尺度なんて間隔尺度でも比尺度でもないので,相関係数なんか計算できないよなあと思われるのでしょう。
しかし,間隔尺度の定義を考えてください。「データの間隔に意味がある」ですね。1と2の間隔と16と17の間隔は同じということです。
特別な場合として,データが2通りの値しか取らない場合,「1と2の間隔」しかないわけで,「じゃあ,だめだろ」といわれるかもしれないけど,「間隔尺度の取りうる値の個数」なんて,どこでも規定されているわけではなく,数学的には「間隔尺度の取りうる値の個数」なんて関係ないのです。二つしか取りうる値がないものも,間隔尺度なんですよ。

蛇足ですが,二値データは,簡単な一次変換でダミー変数になります。

なお,二値データを連続変数と同列に見て統計解析をするのは当たり前のように行われています。ダミー変数をつかっての重回帰分析などは広く行われています。
なお,k ≧ 3値データの場合も k-1 個のダミー変数に展開されて重回帰分析を初めとして多変量回帰に広く使われています。

二値データを数値化するというのには,明らかに意味があるのです。つまり,2つのカテゴリーがあれば,その二つは何らかの意味で大小関係があるのです。「何らかの意味」というのは,そのデータとそのデータを用いた分析の土俵上で決まるものですけど。

黒/白,男/女,高学歴/そうではない, 低収入/高収入,外国人/日本人,ある薬を投与された/されていない/,生産年齢人口/老年人口 ...
倫理的な問題も関与する場合もあるけど,なんだかんだいっても「少なくとも両者には何らかの意味で違いがある」ということは事実なんですから。

「統計学関連なんでもあり」の「全文検索」で「点双列相関係数」を検索すると,いろいろ興味深い記事を読めると思います。

Re: 二値データと連続変数の相関関係について
  投稿者:成瀬 2017/04/23(Sun) 12:46 No. 22328
課長様、青木先生

コメンどうもありがとうございました。
論文を読んでいてもニ値データと連続変数の相関解析の結果を良く目にしていたので、
解析方法に疑問を頂いていました。
2つしか取りうる値がないものも間隔尺度という考えは自分には全くなかったので、
学習になりました。

「点双列相関係数」で検索してさらに勉強させて頂きます。
ありがとうございました。


Cox回帰分析における独立変数の数とイベント数について
  投稿者:BOZ 2017/04/06(Thu) 11:23 No. 22314
連投ですみません。

Cox回帰分析で癌の再発に関連する臨床所見を検討したいと思っています。症例数は172例(うち再発が28例)、臨床所見(独立変数)は107個と非常に多いです。

RでAICによる変数選択(変数増減法)を行い、最終的に4個の独立変数がモデルに選ばれました。しかし以前この掲示板でご紹介いただいた下記の論文によると、イベント数は独立変数の10倍以上ないといけないとされています。ですので私のデータの場合、28/10=2.8なので、独立変数は多く見積もっても3個までとなってしまいます。

https://www.ncbi.nlm.nih.gov/pubmed/8970487?dopt=Abstract&holding=f1000

そこでお聞きしたいのですが、
4個の独立変数を含むモデルは、私のデータの場合ですと不適切でしょうか?
不適切としますと、例えば4個のうち最後に投入された変数1個を除き、独立変数を3個に絞ったモデルとするのは適切でしょうか?

Re: Cox回帰分析における独立変数の数とイベント数について
  投稿者:青木繁伸 2017/04/07(Fri) 22:46 No. 22315
3 はよくて 4 はだめだなんてことがないのは,明らかでしょう。
科学的に考えましょう。
選んだ 3 個も不適切かも知れないと疑うのが先でしょう。

Re: Cox回帰分析における独立変数の数とイベント数について
  投稿者:BOZ 2017/04/07(Fri) 23:13 No. 22316
青木先生、ありがとうございます。
モデルが適切か不適切かはどのようにして確認すればよろしいでしょうか?
モデルが他の集団にも適応可能かどうかバリデーションをすればよろしいですか?

Re: Cox回帰分析における独立変数の数とイベント数について
  投稿者:青木繁伸 2017/04/08(Sat) 19:23 No. 22317
> モデルが適切か不適切かはどのようにして確認

先行研究の知見を総合的に考慮する...

Re: Cox回帰分析における独立変数の数とイベント数について
  投稿者:BOZ 2017/04/08(Sat) 22:53 No. 22318
先行研究の結果と矛盾がなければ、仮にEPV(イベント数÷独立変数の数)≧10を満たさないモデルであっても、不適切とは限らないという事でよろしいでしょうか?

EPV≧10は多変量解析で守るべき指標なのか、あくまで参考程度の指標であるのか、その辺りが統計初心者の私にはよくわからず、悩んでおります。

Re: Cox回帰分析における独立変数の数とイベント数について
  投稿者:BOZ 2017/04/11(Tue) 11:55 No. 22324
EPVが5〜9程度でもいいのではないかという論文もあるようです。
https://www.ncbi.nlm.nih.gov/pubmed/17182981

あまりEPV≧10に固執しないで解析しようと思いました。


ロジスティック回帰分析(AICによる変数選択)の結果について
  投稿者:BOZ 2017/04/06(Thu) 10:36 No. 22313
ロジスティック回帰分析で腫瘍の良性/悪性に関連する画像所見を検討したいと思っています。症例数は70例(良性34例、悪性36例)、画像所見(独立変数)は21個です。

RでAICによる変数選択(変数増減法)を行ったのですが、下記のようにP値が0.05を上回る因子もモデルに入ってしまいます。この場合、有意でなかったFactor_CとFactor_D
を除いたFactor_A+Factor_Bを最終的な良悪性鑑別のモデルとするべきでしょうか?それとも4因子すべて含めるべきでしょうか?
Call:
glm(formula = x[, 1] ~ Factor_A + Factor_B + Factor_C + Factor_D, data = x)

Deviance Residuals:
Min 1Q Median 3Q Max
-0.93851 -0.22333 0.04443 0.19089 0.64343

Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -1.427e+00 9.105e-01 -1.568 0.121822
Factor_A 8.127e-06 2.115e-06 3.843 0.000279 ***
Factor_B 2.625e-01 9.344e-02 2.810 0.006544 **
Factor_C -1.078e+00 6.143e-01 -1.755 0.083930 .
Factor_D -9.557e-01 5.801e-01 -1.648 0.104269
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 0.1024062)

Null deviance: 17.4857 on 69 degrees of freedom
Residual deviance: 6.6564 on 65 degrees of freedom
AIC: 45.947

Number of Fisher Scoring iterations: 2

Re: ロジスティック回帰分析(AICによる変数選択)の結果について
  投稿者:鈴木 康弘 2017/04/09(Sun) 09:13 No. 22319
 Factor_CとFactor_Dを含めないモデルだとAICが45.947より増えてしまう、ということですか?
 だったら4因子とも含めていいんじゃないでしょうか。(あまり自信なし)

Re: ロジスティック回帰分析(AICによる変数選択)の結果について
  投稿者:BOZ 2017/04/09(Sun) 12:36 No. 22320
>>Factor_CとFactor_Dを含めないモデルだとAICが45.947より増えてしまう、ということですか?
そうです。
4因子含めた場合、P値が0.05を上回るのはどう解釈したらよいのでしょうか・・・。

Re: ロジスティック回帰分析(AICによる変数選択)の結果について
  投稿者:青木繁伸 2017/04/09(Sun) 20:36 No. 22321
天はみずからたすくものをたすく
自己責任ってこと
それが研究者が研究者たること
その判断が広く受け入れられるかどうか,研究者としてやっていけるかどうかということ
だれかひとりのお墨付きをえられるかどうということではない
それじゃなきゃ論文一切を提示して私を共著者とする?私がレジェクトすれば,投稿をあきらめる?

Re: ロジスティック回帰分析(AICによる変数選択)の結果について
  投稿者:鈴木 康弘 2017/04/10(Mon) 07:10 No. 22322
AICはモデルの尤度から算出します。
個々の変数が有意かどうかはWald検定だと思います。
尤度比検定とWald検定の結果が一致しない場合、普通は尤度比検定の方を信用するみたいです。

Re: ロジスティック回帰分析(AICによる変数選択)の結果について
  投稿者:BOZ 2017/04/10(Mon) 12:37 No. 22323
青木先生、鈴木様

コメントありがとうございます。
自分の責任にて、4因子すべて含めたモデルで論文投稿してみます。
この度は大変ありがとうございました。


患者さんに4年間介入した結果を分析する方法を教えてください
  投稿者:分析初心者 2017/03/27(Mon) 18:07 No. 22309
初めて投稿いたします。よろしくお願いいたします。
患者さんに介入(入院前と入院中)し、質問紙調査によってデータを収集しました。介入の結果が年々患者さんの意識や態度に変化を及ぼしているのか、いないのかを知りたいと思っております。私たちの介入以外に様々な要因の影響があるとは重々承知しておりますが、検定をしたいと思っております。他に収集している変数は、年齢、性別、診療科、入院病棟などです。
目的変数は、「危険な行動の理解度(5件法)」や「入院中に履いていた履物の種類(3種類、質的データ)」など
意識や行動です。調査は、ある一時点の同じ時期・期間の患者さんを対象に調査しました。しかし、毎年患者さんは異なり、対応のあるデータではありません。介入前、介入1年後、2年後、3年後4年後とそれぞれ400名程度収集しております。このような場合、どのような分析方法が良いのがわかりません。分析方法についてお教えいただけますと幸いです。どうぞよろしくお願いいたします。

Re: 患者さんに4年間介入した結果を分析する方法を教えてください
  投稿者:課長 2017/03/28(Tue) 12:16 No. 22310
「危険な行動の理解度」が介入前後で異なるかどうかなら、マンホイットニーのU検定になるんじゃないでしょうか?経年別に異なるかどうかならクラスカルワリス検定かなと思います。「履物の種類」が介入前後で異なるかどうかなら、カイ二乗検定ですかね。「危険な行動の理解度」や「履物の種類」を従属変数にして、介入や年齢、性別、診療科などのデータを独立変数とした回帰(重回帰や多重ロジスティックなど)も考えられます。いろいろ考えられますが、それぞれにいろいろ問題があるかと思います。

Re: 患者さんに4年間介入した結果を分析する方法を教えてください
  投稿者:分析初心者 2017/03/28(Tue) 17:49 No. 22311
お忙しい中ご助言いただきましてありがとうございます。
>経年別に異なるかどうかならクラスカルワリス検定かなと思います。
経年別にどう変化しているのか、その変化は有意なのかを調べたいと思っております。クラスカルワリス検定は、全く理解しておりませんので、勉強したいと思います。
>いろいろ考えられますが、それぞれにいろいろ問題があるかと思います。
申し訳ありません。その通りだと思っております。
年齢が20歳から95歳までの幅がありましたので、中央値を基準に50才代以下と60才代以上に分けて分析しようと思います。
1年分のデータのみを分析対象としますと目的変数が量的なデータであれば重回帰分析など、「入院中に履いていた履物」のように質的なデータであれば、多重ロジスティック回帰分析ができるということだと理解いたしました。
青木先生の「統計手法の関連図」を順にしていきますと、「クラスカルワリス検定」に行き着きました。
検定の対象が「比率」も「平均値」どちらもありますので、あらゆるものを参考に分析したいと思います。ありがとうございました。

Re: 患者さんに4年間介入した結果を分析する方法を教えてください
  投稿者:課長 2017/03/30(Thu) 01:30 No. 22312
いろいろ大変かと思いますが、がんばってくださいね!


測定回数がばらつくデータの解析方法
  投稿者:大西 2017/03/13(Mon) 09:14 No. 22302
青木先生 
お世話になっております。

現在健診データをもとに、横断研究を予定しています。
2年分の健診データをもとに解析する際に、
2回受診している人と一回しか受診していない人が混在しており、
どのように解析すればよいかわからず困っております。

二回目に受診したデータは除外して、1回目のデータのみ使用する、という方法はもっとも簡単な解決策だとは思いますが、
二回目のデータも活用できるような適切な解析手法がございますでしょうか。

Re: 測定回数がばらつくデータの解析方法
  投稿者:青木繁伸 2017/03/13(Mon) 21:35 No. 22303
> 二回目に受診したデータは除外して、1回目のデータのみ使用する、という方法はもっとも簡単な解決策だとは思いますが、
> 二回目のデータも活用できるような適切な解析手法がございますでしょうか。

一回目しか受信していない人がどれくらいの割合なのかにもよるでしょうが。
そもそも,「二年分の健診データをもとに解析」という意図は,「二回とも受診した人を対象にする」ということなんでしょうか。それとも,漫然と検診して,流れで分析する?

そういうことではないというのならば,必然的に,二回とも受診した人を対象に分析するということになるのでは?

もっとも,二回目しか受診していない人(当人的には1回目なのだが)は対象外になるだろうし,1回目しか受診しなかった人がなぜ2回目を受診しなかったのかというのがバイアスになるかも知れないけど(もう,あんな検診はごめんだ。とか,その後病気になってとか,あるいは死亡したかも知れないし)。

いずれにせよ,事前のデータ設計が不十分だったことの「つけ」をいまさらながらに払わされることになったと言うことですね。

Re: 測定回数がばらつくデータの解析方法
  投稿者:大西 2017/03/14(Tue) 09:07 No. 22304
青木先生

ご返信ありがとうございます。

初回の健診受診時のデータでの解析を基本と考えています。

たとえば、塩分と血圧の関連を見る研究を考えたときに、
ID 健診年度 塩分 血圧
   (g/day) (mmHg)
01 29 15 150
01 30 12 120
02 29 10 110
03 29 9 120
04 29 10 100
04 30 9 100
05 30 12 130
06 30 13   140

29年度を受診した1,2,3,4の29年度分のみを使用するという方法は理解できます。
30年度から初めて受診した5,6番も併せて解析する予定でいました。

ここで、30年度の1,4番のデータは除外されるべきで、解析に加えるべきでないと解釈してよろしいでしょうか。

Re: 測定回数がばらつくデータの解析方法
  投稿者:青木繁伸 2017/03/14(Tue) 14:25 No. 22305
「初回検診時のデータ」が何をあらわしているのでしょうか?
受診者にとって,最初の受診時データなんでしょう?だったら,
> 29年度を受診した1,2,3,4の29年度分のみを使用するという方法は理解できます。
> 30年度から初めて受診した5,6番も併せて解析する予定でいました。
それでよいでしょう。

それ以外のやりかたってあるんですか?

Re: 測定回数がばらつくデータの解析方法
  投稿者:大西 2017/03/15(Wed) 10:53 No. 22307
そうですよね。
大変参考になりました。

ありがとうございました。

Re: 測定回数がばらつくデータの解析方法
  投稿者:tanaka 2017/03/16(Thu) 01:12 No. 22308 HomePage
情報ありがとうございます。


R ftable()関数
  投稿者:明石 2017/03/11(Sat) 19:03 No. 22299
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

先日は、件名「状態遷移を表形式で整理」でご教示いただきまして、誠にありがとうございました。
大変に助かりました、とても良い勉強をさせていただきました。
改めて御礼を申し上げます。

Rプログラムについて、ご教示いただきたいことが出てきました。

何卒どうぞ、よろしくお願いいたします。

ーーー

タイタニックのデータを例に、やりたいことをご説明します。

タイタニックのデータをデータフレーム dat に格納しています。

> dim(dat)
[1] 2201 4

> colnames(dat)
[1] "Class" "Age" "Sex" "Alive"

> ftable(dat)

Alive 死亡 生還
Class Age Sex
1等 子ども 女性 0 1
男性 0 5
大人 女性 4 140
男性 118 57
2等 子ども 女性 0 13
男性 0 11
大人 女性 13 80
男性 154 14
3等 子ども 女性 17 14
男性 35 13
大人 女性 89 76
男性 387 75
乗務員 子ども 女性 0 0
男性 0 0
大人 女性 3 20
男性 670 192

左側のカラムから順にソートされて整形されます。

もし、以下のように、右側からソートされるような整形をしたいと思った場合に、
どのようにすればよろしいのでしょうか?

Alive 死亡 生還
Class Age Sex
1等 大人 男性 118 57
2等 大人 男性 154 14
3等 大人 男性 387 75
乗務員 大人 男性 670 192
1等 子ども 男性 0 5
2等 子ども 男性 0 11
3等 子ども 男性 35 13
乗務員 子ども 男性 0 0
1等 大人 女性 4 140
2等 大人 女性 13 80
3等 大人 女性 89 76
乗務員 大人 女性 3 20
1等 子ども 女性 0 1
2等 子ども 女性 0 13
3等 子ども 女性 17 14
乗務員 子ども 女性 0 0

色々と考えていますが、できません。
ご教示いただければ大変に助かります。

何卒どうぞ、よろしくお願いいたします。

Re: R ftable()関数
  投稿者:青木繁伸 2017/03/11(Sat) 23:27 No. 22300
引数指定で処理できるようではなく,また,既存の関数もないので,自分で書くしかないかな??

結果のヘッダ部分が ftable のものと違うが,以下のような感じで
> a = ftable(dat)
> b = data.frame(a)
> # ftable の結果をデータフレームにして,右の列から並べ替える
> # 要素が factor のときは,decreasing を要素ごとに指定することができないので,as.integer で整数コードの方を用いる
> # method はデフォルトではなく "radix" を指定すること
> # 並べ替えのための以下の 1 行は,それぞれのデータフレームにより適切に記述のこと
> result = b[order(as.integer(b$Sex), as.integer(b$Age), as.integer(b$Class),
+ decreasing=c(TRUE, TRUE, FALSE),
+ method="radix"),]
> result
Class Age Sex Alive Freq
13 1等 大人 男性 死亡 118
29 1等 大人 男性 生還 57
14 2等 大人 男性 死亡 154
30 2等 大人 男性 生還 14
15 3等 大人 男性 死亡 387
31 3等 大人 男性 生還 75
16 乗務員 大人 男性 死亡 670
32 乗務員 大人 男性 生還 192
9 1等 子ども 男性 死亡 0
25 1等 子ども 男性 生還 5
10 2等 子ども 男性 死亡 0
26 2等 子ども 男性 生還 11
11 3等 子ども 男性 死亡 35
27 3等 子ども 男性 生還 13
12 乗務員 子ども 男性 死亡 0
28 乗務員 子ども 男性 生還 0
5 1等 大人 女性 死亡 4
21 1等 大人 女性 生還 140
6 2等 大人 女性 死亡 13
22 2等 大人 女性 生還 80
7 3等 大人 女性 死亡 89
23 3等 大人 女性 生還 76
8 乗務員 大人 女性 死亡 3
24 乗務員 大人 女性 生還 20
1 1等 子ども 女性 死亡 0
17 1等 子ども 女性 生還 1
2 2等 子ども 女性 死亡 0
18 2等 子ども 女性 生還 13
3 3等 子ども 女性 死亡 17
19 3等 子ども 女性 生還 14
4 乗務員 子ども 女性 死亡 0
20 乗務員 子ども 女性 生還 0
>
> # これ以降は,データフレームの列名等に左右されないように汎用的に記述する
> name = levels(result[,ncol(result)-1]) # この時点では,result の右から2列目が目的の列名(例では "Alive")
> # name はその水準名(例では "死亡","生還"),これを後で使う
> lname = length(name) # lname はその長さ
> Freq = matrix(result[,"Freq"], byrow=TRUE, ncol=lname) # 集計数がある列から数値を取り出し,水準数を列とする行列にする
> result = cbind(result[seq(1, nrow(result), by=lname), 1:(ncol(result)-2)], Freq) # 行名と集計数の行列をバインド
> names(result)[ncol(result)-lname:1+1] = name # 集計数の列名を付ける
> print(result, row.names=FALSE)
Class Age Sex 死亡 生還
1等 大人 男性 118 57
2等 大人 男性 154 14
3等 大人 男性 387 75
乗務員 大人 男性 670 192
1等 子ども 男性 0 5
2等 子ども 男性 0 11
3等 子ども 男性 35 13
乗務員 子ども 男性 0 0
1等 大人 女性 4 140
2等 大人 女性 13 80
3等 大人 女性 89 76
乗務員 大人 女性 3 20
1等 子ども 女性 0 1
2等 子ども 女性 0 13
3等 子ども 女性 17 14
乗務員 子ども 女性 0 0

しかしまあ,row.vars と col.vars を使って,以下のようにした方が簡単でしょう?
Class を一番左に置きたいという理由がないならば
男性と女性,大人と子どもの前後関係を変えるのは,factor でまえもって処理しておけばよいだけ,だし。
> ftable(dat, row.vars=c("Sex", "Age", "Class"), col.vars="Alive")
Alive 死亡 生還
Sex Age Class
女性 子ども 1等 0 1
2等 0 13
3等 17 14
乗務員 0 0
大人 1等 4 140
2等 13 80
3等 89 76
乗務員 3 20
男性 子ども 1等 0 5
2等 0 11
3等 35 13
乗務員 0 0
大人 1等 118 57
2等 154 14
3等 387 75
乗務員 670 192

Re: R ftable()関数
  投稿者:明石 2017/03/12(Sun) 11:45 No. 22301
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

ご丁寧に、色々な方法をご教示くださいまして、誠にありがとうございます。
とても勉強になります。

最初、ftable の結果をデータフレームにして並べ替えることを検討しましたが、
うまくいきませんでした。

青木先生がお示しをしてくださいましたプログラムを見て、理解できました。

心から御礼を申し上げます。
ありがとうございました。


サンプルサイズ
  投稿者:Bob 2017/03/10(Fri) 21:01 No. 22296
青木先生
お世話になっております。
わたくし、統計の初歩がわかってない者でございます。
二群でトレーニング日数の違いがあるのかを調べる臨床研究を行うため、過去の報告に基づき、サンプルサイズの計算をしようとしていますが、その方法がわからず閉口しています。

過去の類似の報告で、患者さんのトレーニング日数を二群で比較したものがございます。
その結果に基づき、サンプルサイズを算出しようと思っています。

過去の報告:
A群のトレーニング日数:5.4日±2.7日
B群のトレーニング日数:3.6日±1.4日

そこで以下についてお教え頂けないでしょうか。
t検定を行う前提で、
サンプルサイズを計算するとき
予め、両群の日数の差と標準偏差が必要ということを理解しています。
この時の標準偏差を求める方法をお教え頂けないでしょうか。
恐らく、±2.7日と±1.4日の併合標準偏差を求めないといけないのではないかと察していますが、種々しらべてもその方法がわかりませんでした。
わたくしの理解が正しければ併合標準偏差の求め方をお教え頂ければ幸いでござます。

または、適切な標準偏差の求め方をお教え頂ければ幸いでござます。

基本的な質問で恐縮でございますが、どうぞ宜しくお願い申し上げます。

Re: サンプルサイズ
  投稿者:青木繁伸 2017/03/10(Fri) 22:06 No. 22297
http://aoki2.si.gunma-u.ac.jp/lecture/Univariate/combine.html
をみれば,一般的な場合の解法が書いてあります。
二群の場合はより簡単に
http://aoki2.si.gunma-u.ac.jp/lecture/Average/t-test.html
の Ue を計算する式で示してあります。

しかし,サンプルサイズの計算における標準偏差は,理論的な基礎はあるものの数値的にはストリクトなものととらえる必要はないしそのようにすべきものでもないでしょう。なんといっても,その計算値はある特定の場合の数値に依存しているのだから,その基礎となった数値の振れ幅(誤差)も考慮しなければならないと言うことでしょう。
例えば,上記のように計算された数値を踏まえて,実際的にはその前後の何個かの値を用いて,最も厳しい場合,楽観的な場合というような推計値を計算し,さらにそれらを総合的に判断するというようなことが必要なのではないでしょうか??

Re: サンプルサイズ
  投稿者:Bob 2017/03/10(Fri) 22:38 No. 22298
青木先生
ご多忙の折、早々とご回答ありがとうございました。
とても勉強になりました。
とても感謝しております。
お礼まで


R 状態遷移を表形式で整理(2)
  投稿者:明石 2017/03/05(Sun) 17:31 No. 22293
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

昨日〜本日と、件名「状態遷移を表形式で整理」について、
有難いご教示をいただきまして、誠にありがとうございました。
御礼を申し上げます。

その継続について検討を進めています。

私の方でもプログラムを作成していますが、苦慮しております。
またご教示いただれば、大変に助かります。

何卒どうぞ、よろしくお願いいたします。

ーーーーーーー

今回は、Webページのコンテンツの遷移を、
・縦方向に、遷移前のコンテンツ名
・横方向に、遷移後のコンテンツ名
・セルに、出現度数
からなる表形式で整理したいと思います。

今回の例示では、コンテンツの状態遷移が17行です、

n05 > n07 > n01 > n05
n05 > n05 > n05
n02 > n05
n07 > n02 > n05 > n02 > n05
n05
n03 > n05
n07 > n07 > n02
n02 > n02 > n04 > n06 > n01
n02 > n01 > n02 > n01 > n04 > n02
n06 > n02 > n02 > n06 > n01 > n02 > n06 > n01
n07 > n06 > n06
n04 > n03
n07 > n05 > n01 > n01
n02 > n04 > n04 > n02
n05 > n06 > n01 > n02
n05 > n02 > n02
n02 > n07 > n01 > n01 > n01

これを読み取り、画像ファイルにお示しをします表を出力したいと思います。

以下の手順でプログラムを作成するつもりです。
誠に恥ずかしいのですが、リストのハンドリングが不得手ですので、最初の(1)で足踏み状態です。

(1)
各行の状態遷移を読み取り、遷移の前後の対を出力する。

例えば、
"n05 > n07 > n01 > n05" では、
"n05_n07" , "n07_n01" , "n01_n05" のように、
遷移前のコンテンツ、遷移後のコンテンツを、"_"で連結した文字列を出力します。

ただし、以下の2つの場合は、無視します。。
・ "n05 > n05" のように、同じコンテンツが連続する場合
・ "n05" のように、1行でコンテンツが1つの場合


(2)
青木先生からご教示いただきましたと同じ方法で、表を出力します。

縦方向は、セッション
横方向は、遷移の前後の対の文字列

(3)
上記の表で、行和が 0 の行を削除する。

つまり、以下の場合です。
・1行に1つのコンテンツしかない場合
・1行に複数のコンテンツがあるが、すべて同じ種類の場合

(4)
上記の表で、列和をとる。

(5)
上記の表から、列名と列和の2列からなるデータフレームを作成する。

2列の変数名です。
mae_ato, freq

(6)
上記のデータフレームで、
変数 mae_ato を strsplit() 関数で分割した、3列を作成する。

3列の変数名です。
mae, ato, freq

(7)
上記のデータフレームで、xtab()関数を用いて、クロス表を作成する。

縦方向は、mae
横方向は、ato
セルは、freq

これが最終的に求める表です。

青木先生は、研ぎ澄まされた本質的でスマートなプログラムを作成されますので、
ご教示いただけましたら、大変に良い勉強となります。

お手数をお掛けいたしますが、ご教示を頂戴できましたら助かります。

何卒どうぞ、よろしくお願いいたします。


Re: R 状態遷移を表形式で整理(2)
  投稿者:青木繁伸 2017/03/05(Sun) 21:24 No. 22294
以下のように(前のプログラム方針とちょっと違うところがあるが,まあそれはそれ)
func <- function(s) {
s <- strsplit(s, " > ")
nameOfContents <- names(table(unlist(s)))
numberOfContents <- length(nameOfContents)
tbl <- matrix(0, numberOfContents, numberOfContents)
colnames(tbl) <- rownames(tbl) = nameOfContents
for (S in s) {
a <- unlist(S)
n <- length(a)
for (i in seq_len(n - 1)) {
if (a[i] != a[i + 1]) {
tbl[a[i], a[i + 1]] <- tbl[a[i], a[i + 1]] + 1
}
}
}
tbl
}

func(c("n05 > n07 > n01 > n05",
"n05 > n05 > n05",
"n02 > n05",
"n07 > n02 > n05 > n02 > n05",
"n05",
"n03 > n05",
"n07 > n07 > n02",
"n02 > n02 > n04 > n06 > n01",
"n02 > n01 > n02 > n01 > n04 > n02",
"n06 > n02 > n02 > n06 > n01 > n02 > n06 > n01",
"n07 > n06 > n06",
"n04 > n03",
"n07 > n05 > n01 > n01",
"n02 > n04 > n04 > n02",
"n05 > n06 > n01 > n02",
"n05 > n02 > n02",
"n02 > n07 > n01 > n01 > n01"))
実行結果
    n01 n02 n03 n04 n05 n06 n07
n01 0 3 0 1 1 0 0
n02 2 0 0 2 3 2 1
n03 0 0 0 0 1 0 0
n04 0 2 1 0 0 1 0
n05 1 2 0 0 0 1 1
n06 4 1 0 0 0 0 0
n07 2 2 0 0 1 1 0

Re: R 状態遷移を表形式で整理(2)
  投稿者:明石 2017/03/06(Mon) 09:29 No. 22295
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

今回も、研ぎ澄まされたプログラムをご教示いただきまして、誠にありがとうございました。

確認いたしました。

ご教示くださいました2種類のプログラム、しっかりと勉強させていただきます。

心から御礼を申し上げます。


R 状態遷移を表形式で整理
  投稿者:明石 2017/03/04(Sat) 18:46 No. 22288
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

Rプログラムについて、ご教示いただきたいことが出てきました。

何卒どうぞ、よろしくお願いいたします。

ーーー

Webページのコンテンツの遷移を、表形式で整理する、Rプログラムを作成したいと思います。

簡単な例でお示しをします。

n05 > n07 > n01 > n05
n07 > n07 > n07 > n07
n02 > n05
n07 > n02 > n05 > n02 > n05
n05
n07 > n07 > n02
n07 > n05 > n01 > n01
n02 > n07 > n01 > n01 > n01
・・・・・・・

各行は、1つのセッションです。

n01 n02 n03,,,n07 の文字列は、コンテンツの名前です。

例えば、1行目は、n05 → n07 → n01 → n05 の4つのコンテンツを閲覧したという意味です。

2行目のように、同じコンテンツだけという場合もありえます。

この状態遷移を、画像でお示しします表の形で整理したいと思います。

縦方向は、セッション
横方向は、コンテンツの種類

セルには、各セッションで閲覧したコンテンツ種類ごとの回数を記入します。

実際のデータでは、データファイルごとに、コンテンツの名称も個数も変わります。

固定ではないところが、プログラム作成の難しさの要因ともなっています。

ご教示いただけましたら、大変に助かります。
何卒どうぞ、よろしくお願いいたします。


Re: R 状態遷移を表形式で整理
  投稿者:青木繁伸 2017/03/04(Sat) 19:48 No. 22289
コンテンツ名に特定の規則がない場合にも対応できるように
臆せず for を使う
func <- function(s) {
tbl <- NULL
s <- strsplit(s, " > ")
nameOfContents <-(table(unlist(s)))
tbl <- matrix(0, length(s), length(nameOfContents))
colnames(tbl) <- nameOfContents
for (i in seq_along(s)) {
for (j in seq_along(s[[i]])) {
tbl[i, s[[i]][j]] <- tbl[i, s[[i]][j]]+1
}
}
tbl
}

func(c("n05 > n07 > n01 > n05",
"n07 > n07 > n07 > n07",
"n02 > n05",
"n07 > n02 > n05 > n02 > n05",
"n05",
"n07 > n07 > n02",
"n07 > n05 > n01 > n01",
"n02 > n07 > n01 > n01 > n01"))
実行結果
     n01 n02 n05 n07
[1,] 1 0 2 1
[2,] 0 0 0 4
[3,] 0 1 1 0
[4,] 0 2 2 1
[5,] 0 0 1 0
[6,] 0 1 0 2
[7,] 2 0 1 1
[8,] 3 1 0 1


Re: R 状態遷移を表形式で整理
  投稿者:明石 2017/03/05(Sun) 08:55 No. 22290
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

ご教示をいただきまして、誠にありがとうございます。

追試しながら、勉強をさせていただいております。
内容については、理解しつつあります。

しかしながら、掲示板のプログラムをコピペさせていただき、実行すると、
以下のようなエラー表示があります。

誠にお手数ですが、ご確認をいただけると大変に助かります。

どうぞ、よろしくお願いいたします。

ご教示をいただきまして、誠にありがとうございます。

> func <- function(s) {
+ tbl <- NULL
+ s <- strsplit(s, " > ")
+ nameOfContents <-(table(unlist(s)))
+ tbl <- matrix(0, length(s), length(nameOfContents))
+ colnames(tbl) <- nameOfContents
+ for (i in seq_along(s)) {
+ for (j in seq_along(s[[i]])) {
+ tbl[i, s[[i]][j]] <- tbl[i, s[[i]][j]]+1
+ }
+ }
+ tbl
+ }
>
> func(c("n05 > n07 > n01 > n05",
+ "n07 > n07 > n07 > n07",
+ "n02 > n05",
+ "n07 > n02 > n05 > n02 > n05",
+ "n05",
+ "n07 > n07 > n02",
+ "n07 > n05 > n01 > n01",
+ "n02 > n07 > n01 > n01 > n01"))

tbl[i, s[[i]][j]] でエラー: 添え字が許される範囲外です
>

Re: R 状態遷移を表形式で整理
  投稿者:青木繁伸 2017/03/05(Sun) 09:39 No. 22291
4 行目がおかしかったようです(関数名 names が落ちてしまった)
nameOfContents <- names(table(unlist(s)))

Re: R 状態遷移を表形式で整理
  投稿者:明石 2017/03/05(Sun) 10:13 No. 22292
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

ご教示をいただきまして、誠にありがとうございます。

何度もお手を煩わせてしまい、申し訳ございませんでした。

無事に確認できました。
プログラムの内容も、理解できました。
ここで使われているノウハウはとても貴重であり、色々と役立ちます!

お示しくださいましたプログラムは、
状態遷移を記載した外部ファイル"Transition.txt"を、
scan()関数で読み込むように、汎用化させていただきました。

ご教示をいただきまして、誠にありがとうございます。
心から御礼を申し上げます。


fisherの正確検定について
  投稿者:医学系の大学院生 2017/02/26(Sun) 14:21 No. 22285
青木様

いつも拝見させていただいております。
医学系の研究における統計手法について質問がございます。

私は現在、ある上皮癌を浸潤群と非浸潤群の2群間について治療法や診断法についての検討を行っています。
この上皮癌には複数の組織型(17種類)が含まれます。
組織型 浸潤群 非浸潤群
 
A型 1 0
B型 0 1
C型 2 4
D型 2 3
E型 1 9
F型 3 0
G型 1 0
H型 0 1
I型 17 31
J型 0 3
K型 0 7
L型 2 1
M型 15 67
N型 3 0
O型 4 4
P型 3 2
Q型 3 6
2群間における組織型の症例数の割合に差があるかどうかを調べたいので、2x17のfisherの正確検定を行いました。統計ソフトはJMPを用いて行ってみたのですが、計算量が多いためか、途中でフリーズしてしまいます。また、先生のホームページhttp://aoki2.si.gunma-u.ac.jp/exact/fisher/fisher.cgiを用いた場合も、計算中です…の表示のままです。

検定はfisherの正確検定でよいのか、また、上記の解決方法などございましたらご助言いただければ幸いです。

Re: fisherの正確検定について
  投稿者:青木繁伸 2017/02/26(Sun) 20:33 No. 22286
R を使ってください
> b
浸潤群 非浸潤群
A型 1 0
B型 0 1
C型 2 4
D型 2 3
E型 1 9
F型 3 0
G型 1 0
H型 0 1
I型 17 31
J型 0 3
K型 0 7
L型 2 1
M型 15 67
N型 3 0
O型 4 4
P型 3 2
Q型 3 6

> fisher.test(b, workspace=2000000)

Fisher's Exact Test for Count Data

data: b
p-value = 0.0003557
alternative hypothesis: two.sided

Re: fisherの正確検定について
  投稿者:医学系の大学院生 2017/02/26(Sun) 21:14 No. 22287
青木様

早速のアドバイスありがとうございます。
Rについては使用したことがなかったので、再度勉強したのちに統計をとってみようと思います。
ありがとうございました。


距離行列からの特定データ抽出
  投稿者:shimp 2017/02/24(Fri) 15:51 No. 22282
Rを用いて得た距離行列から特定の条件を満足する組合せを抽出したいと考えています。

例えば、
(USA.dist <- dist(USArrests)^2)
で得られる距離行列から、値10,000を超える組合せのみを抽出したいと考えています。
抽出するためには、データフレームへの変換が必要だと考えて、
(USA.dist <- as.data.frame(as.matrix(USA.dist)))
としたうえで、
(subset(USA.dist, USA.dist > 10000))
としてみたのですが、NAで埋まった行名も失われた出力となってしまいます。

subsetの書き方がまずい様な気がするのですが、そもそもas.data.frame(as.matrix())といった使い方がいかにも力づくだと感じています。
どうかアドバイスを賜りたく何卒宜しくお願い致します。

Re: 距離行列からの特定データ抽出
  投稿者:青木繁伸 2017/02/25(Sat) 09:48 No. 22283
「特定データの抽出」とは,どのようなことを意味するのでしょうか?
> (x <- c(a=2, b=4, c=1, d=6, e=7))
a b c d e
2 4 1 6 7

> (y <- dist(x)^2) # class(y) は "dist" 三角行列
a b c d
b 4
c 1 9
d 16 4 25
e 25 9 36 1

> (z <- as.matrix(y)) # class(z) は "matrix" 正方行列
a b c d e
a 0 4 1 16 25
b 4 0 9 4 9
c 1 9 0 25 36
d 16 4 25 0 1
e 25 9 36 1 0

> # ある値より大きい要素の添え字
> which(z > 20) # デフォルトでは,arr.ind=FALSE になっており,結果はベクトルでの index
[1] 5 14 15 18 21 23

> (index <- which(z > 20, arr.ind=TRUE)) # arr.ind= TRUE で,行列での index
row col
e 5 1
d 4 3
e 5 3
c 3 4
a 1 5
c 3 5

> # 実際のその要素
> z[index]
[1] 25 25 36 25 25 36 # ベクトルからの抽出

> # 行列のまま抽出しようという試みは挫折する(これでもいいということか?)
> z[index, , drop=FALSE]
a b c d e
e 25 9 36 1 0
d 16 4 25 0 1
e 25 9 36 1 0
c 1 9 0 25 36
a 0 4 1 16 25
c 1 9 0 25 36
a 0 4 1 16 25
c 1 9 0 25 36
c 1 9 0 25 36
d 16 4 25 0 1
e 25 9 36 1 0
e 25 9 36 1 0

> # ある値より大きい要素の添え字を行列名に変換
> matrix(colnames(z)[t(index)], byrow=TRUE, ncol=2)
[,1] [,2]
[1,] "e" "a"
[2,] "d" "c"
[3,] "e" "c"
[4,] "c" "d"
[5,] "a" "e"
[6,] "c" "e"

> # これではないのでしょうけど
> u <- z
> u[u < 20] <- 0
> u
a b c d e
a 0 0 0 0 25
b 0 0 0 0 0
c 0 0 0 25 36
d 0 0 25 0 0
e 25 0 36 0 0

> # こっちとか
> v <- z
> b <- apply(v, 1, function(a) any(a > 20))
> v[b, b]
a c d e
a 0 1 16 25
c 1 0 25 36
d 16 25 0 1
e 25 36 1 0

Re: 距離行列からの特定データ抽出
  投稿者:shimp 2017/02/25(Sat) 15:20 No. 22284
有難う御座います。説明が不適切・不足していたことをお詫びします。

やりたかったことは、まさしく最後のケースです。
また最初の質問で" > 10000"としていたのは誤りで、やりたいことはご近所探しであったため不等号の向きが逆でした。重ねて申し訳ありません。

5,000件くらいの膨大なデータの中からご近所を取り出して、距離の数字を確認しながらデンドログラムを眺めてみたいと思ったのです。

par(mfrow = c(1, 2))

x <- c(a=1, b=2, c=5, d=8, e=11, f=14, g=17, h=20, i=23, j=26, k=29, l=32, m=35, n=38, o=41, p=42, q=43)
(y <- dist(x)^2)
plot(hclust(y),hang=-1)

v <- as.matrix(y)
b <- apply(v, 1, function(a) any(a < 3 & a > 0))
(w <- as.dist(v[b, b]))
plot(hclust(w),hang=-1)

週明けに実際のデータに当ててみるのが楽しみでなりません。
有難う御座いました。


AIC による,ヒストグラム(度数分布表)の最適階級分割の探索
  投稿者:明石 2017/02/23(Thu) 11:28 No. 22277
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

Rプログラムについて、ご教示いただきたいことが出てきました。

何卒どうぞ、よろしくお願いいたします。

ーーー

http://aoki2.si.gunma-u.ac.jp/R/AIC-Histogram.html
AIC による,ヒストグラム(度数分布表)の最適階級分割の探索

を利用させていただいております。
大変に重宝しております。御礼を申し上げます。

情報量基準では、AIC の他に BIC があります。

基準が違うので、それぞれの基準での最適化が図られており、
その違いを確認したいと思っております。

私が扱うデータには、年齢、年収などの項目があり、外れ値がありますので、
基準の違いによる最適化の違いに強い興味があります。

しかしながら、AIC、BIC の公式を本で見たことはありますが、
Rでプログラムを書けるまでの知識はありません。

誠に厚かましいお願いでございますが、

AIC による,ヒストグラム(度数分布表)の最適階級分割の探索

BIC による,ヒストグラム(度数分布表)の最適階級分割の探索
のRプログラムをご教示いただければ、大変に助かります。

どうぞよろしくお願いいたします。
失礼いたします。

Re: AIC による,ヒストグラム(度数分布表)の最適階級分割の探索
  投稿者:青木繁伸 2017/02/23(Thu) 18:01 No. 22278
R の関数 AIC のオンラインヘルプには,
Generic function calculating Akaike's ‘An Information Criterion’ for one or several fitted model objects for which a log-likelihood value can be obtained, according to the formula -2*log-likelihood + k*npar, where npar represents the number of parameters in the fitted model, and k = 2 for the usual AIC, or k = log(n) (n being the number of observations) for the so-called BIC or SBC (Schwarz's Bayesian criterion).
と書いてある。
要するに -2*log-likelihood + k*npar において,AIC では k = 2,BIC では log(サンプルサイズ)
なので,AIC.Histogram の
AIC <- -2*(logNZ(sum1, c1*N)+temp+logNZ(sum2, c2*N))+2*cc1c2r1

AIC <- -2*(logNZ(sum1, c1*N)+temp+logNZ(sum2, c2*N))+log(N)*cc1c2r1
として,AIC という文字列を BIC に全置換すればよいと思われ...

【御礼】BIC による,ヒストグラム(度数分布表)の最適階級分割の探索
  投稿者:明石 2017/02/23(Thu) 21:13 No. 22279
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

ご教示をいただきまして、誠にありがとうございます。

AIC、BICでは、基準が違うので、
それぞれの基準での最適化が図られており、
その違いを確認したいと思っております。

ありがとうございました。
心から御礼を申し上げます。

Re: AIC による,ヒストグラム(度数分布表)の最適階級分割の探索
  投稿者:青木繁伸 2017/02/23(Thu) 22:35 No. 22280
log(N) は定数なので,この場合は AIC と BIC は(サンプルサイズ N) が同じなので,平行移動なだけでは?つまり,モデル選択に差はない(??ちょっとおかしいか??)

Re: AIC による,ヒストグラム(度数分布表)の最適階級分割の探索
  投稿者:明石 2017/02/24(Fri) 10:32 No. 22281
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

最初は、AIC , BIC では結果が異なるだろうと思っていました。

手持ちのデータ(約5000件、外れ値を含む)で、色々とやってみましたが、
ヒストグラムの bin の数、bin の区間幅も、一致しました。

カンニングをした後付けですが、
この題材では、青木先生がお書きになられた通りの理解です。

良い勉強をさせていただきました。
改めて御礼を申し上げます。


一般化線形モデルにおける多重比較法
  投稿者:太郎 2017/02/07(Tue) 08:52 No. 22273
以下は引用データですが、TRTが3種類の薬剤処理、REPが反復、Disは死滅株数、Healthは健全株数の1因子乱塊法実験です。
一般化線形モデルを用いて以下のように、Tukey法で多重比較を行いました。これをTukeyではなく、普通の線形モデルにおける分散分析後のLSD法のような方法で比較したいのですが、どなたかご存知の方がおられましたらご教授いただけないでしょうか。

TRT<-c("a","a","a","b","b","b","c","c","c")
REP<-c("1","2","3","1","2","3","1","2","3")
Dis<-c(34,30,28,14,13,24,21,28,25)
Health<-c(12,12,9,28,36,22,24,16,15)
Dataset<-data.frame(TRT,REP,Dis,Health)
library(car)
GLM <- glm(cbind(Dis,Health) ~ TRT + REP , family=binomial(logit), data=Dataset)
Anova(GLM, type="II")
summary(GLM)
library(multcomp)
tukey1<- glht(GLM,linfct=mcp(TRT = "Tukey"))
summary(tukey1)

Re: 一般化線形モデルにおける多重比較法
  投稿者:太郎 2017/02/09(Thu) 09:02 No. 22276
 レスがつかないので、やはり難しいのですね。
 あまりスマートな方法ではないのですが、2つの処理ごとに以下のような分析を行うことにしました。
 普通の線形モデルで、2つの処理ごとにt検定を繰り返し適用するような感じですね。

Dataset<-data.frame(
TRT<-c("a","a","a","c","c","c"),
REP<-c("1","2","3","1","2","3"),
Dis<-c(34,30,28,21,28,25),
Health<-c(12,12,9,24,16,15)
)

library(car)
GLM<-glm(cbind(Dis,Health)~TRT+REP,family=binomial(logit),data=Dataset)
Anova(GLM,type="II")


ポアソン分布への適合度検定におけるカテゴリー併合時の期待度の値
  投稿者:いわさき 2017/02/05(Sun) 21:44 No. 22270
いつも拝見させて頂いております。

適合度の検定 --ポアソン分布への適合度の検定
http://aoki2.si.gunma-u.ac.jp/lecture/GoodnessOfFitness/poissondist.html

にて、
4.期待値が 1 以下のカテゴリーを併合する
とされておられます。通常(他のWEB情報など) 5 以下が多いのですが、
どのようなご判断から、1とすることができるのでしょうか?

現在同様の検定、n=370、平均0.1(推定)、にてカテゴリー値 5 以下を併合してしまうと
自由度が0となり、カイ2乗検定ができなくってしまい悩んでおります。

どうぞよろしくお願い致します。

Re: ポアソン分布への適合度検定におけるカテゴリー併合時の期待度の値
  投稿者:青木繁伸 2017/02/06(Mon) 04:38 No. 22271
世の中には,不正確なものや,自分で確認することなく誤った情報の孫引きがはびこっているので,注意が必要。

Cochran, W. G.: Some methods for strengthening the common χ2 tests. Biometrics, 10, 417-451, 1954
この中で,クロス集計表の独立性の検定の際の最小期待値についてちゃんと書かれているのは有名
http://aoki2.si.gunma-u.ac.jp/lecture/Cross/warning.html
だが,適合度の検定についても書かれている。
つまり,独立性の検定の場合も,最小期待値は「1以上」


Re: ポアソン分布への適合度検定におけるカテゴリー併合時の期待度の値
  投稿者:いわさき 2017/02/06(Mon) 23:22 No. 22272
青木先生

ご指導ありがとうございました。
また調査不足にて大変失礼致しました。ご記述の論文を入手致します。
それを読み終えてから、追加のご質問すべきとは存じますが、
もし可能でしたら、もう一点ご一言頂けますと幸甚です。

次のような、ポアソン分布への適合度を検定において、
お教えの英文から、カテゴリー併合は期待値aでOKのように読みましたが、
さらに、期待値bまで行うべきでしょうか?

再々の初歩的質問にて大変恐縮です。何卒よろしくお願い申し上げます。

==ここから、適合度計算のためのカテゴリー併合の流れ
実測値
階級 0 1 2 3
実測値  340 31 5 0

平均(推定)≒0.11の期待値
階級 0 1 2 3
期待値 336.8 37.1 2.0 0.1

期待値が1未満の桝目を併合
階級 0 1 2
期待値a 336.8 37.1 2.1

さらに、
期待値aに5未満の桝目が全体の桝目の数の20%以上あるので
階級 0 1
期待値b 336.8 39.2
実測値  340 36
しかし自由度0となり検定不能。

===ここまで

Re: ポアソン分布への適合度検定におけるカテゴリー併合時の期待度の値
  投稿者:青木繁伸 2017/02/07(Tue) 09:57 No. 22274
期待値5以下のセルが 20% 以下というのは,k×l 分割表のときの話

Re: ポアソン分布への適合度検定におけるカテゴリー併合時の期待度の値
  投稿者:いわさき 2017/02/07(Tue) 20:53 No. 22275
ご指導ありがとうございました。


因子分析 主因子法と最尤法
  投稿者:kashio 2017/01/31(Tue) 12:28 No. 22266
こんにちは

いつもお世話になります。現在は因子分析バリマックス回転は下位尺度の相関関係が0と仮定できなけラバ用いない。ということで、よくやられている最尤法プロマックス回転でやってみたところ、二重負荷が半数以上に見られたので主因子法プロマックス回転をしたら、全く2重不可はかからず、いい結果だと思いました。そこで、最尤法と主因子法何が違くてこんなことになるのか、そして、今はやりの最尤法ではなく主因子法で問題ないのかご助言ください。

Re: 因子分析 主因子法と最尤法
  投稿者:青木繁伸 2017/01/31(Tue) 13:22 No. 22267
最尤法と主因子法の違いの原因とか,「今はやり」のものでなくてよいのかには答えませんが,R の因子分析では psych ライブラリの fa 関数が有名ですが,因子抽出法としては,以下の 7 種類もあり,しかも minres がデフォルトです(最尤法 ml ではありません)。
factoring method
 fm="minres" will do a minimum residual (OLS),
 fm="wls" will do a weighted least squares (WLS) solution,
 fm="gls" does a generalized weighted least squares (GLS),
 fm="pa" will do the principal factor solution,
 fm="ml" will do a maximum likelihood factor analysis.
 fm="minchi" will minimize the sample size weighted chi square when treating pairwise correlations with different number of subjects per pair.
 fm ="minrank" will do a minimum rank factor analysis.
ちなみに,回転法もたくさんあります。斜交回転でも 8 種類。デフォルトは oblimin です。
orthogonal rotations
 "none", "varimax", "quartimax", "bentlerT",
 "equamax", "varimin", "geominT", "bifactor"
oblique transformations
 "Promax", "promax", "oblimin", "simplimax",
 "bentlerQ, "geominQ", "biquartimin", "cluster"

The default is to do a oblimin transformation, although versions prior to 2009 defaulted to varimax.
SPSS seems to do a Kaiser normalization before doing Promax, this is done here by the call to "promax" which does the normalization before calling Promax in GPArotation.

「最尤法でプロマックス」に限る必要はないし,むしろそれに限らない方がよいのではということかな?

Re: 因子分析 主因子法と最尤法
  投稿者:kashio 2017/01/31(Tue) 18:06 No. 22268
ありがとうございました。デファルトでもやってみます。


ホスマーレメショウ検定のp値について
  投稿者:浅井 2017/01/29(Sun) 09:24 No. 22263
初歩的なことで申しわけありませんが、ご質問させていただきます。
ロジスティック回帰分析のホスマーレメショウ検定のp値についてです。
SPSSにて解析し、ホスマーレメショウ検定以外のp値は有意に問題なく出ていたのですが、ホスマーレメショウ検定のみp値が「.」と表示されていました。これはどのように解釈すればよろしいのでしょうか?ちなみにnは30vs15にて行いました。

Re: ホスマーレメショウ検定のp値について
  投稿者:青木繁伸 2017/01/29(Sun) 10:36 No. 22264
データ数が少なすぎるんじゃないですか?
自分で手計算してみると,その過程で「これじゃ検定できないよなあ」と気づくかも知れませんよ。
せめて,出力結果をそのまま見せてください。

Re: ホスマーレメショウ検定のp値について
  投稿者:浅井 2017/01/29(Sun) 10:58 No. 22265
早速のご回答ありがとうございます。
今、資料が手元にないので後日送信させていただきますが、今後データ数が増える予定がありません。
このような結果では、発表は難しいですよね…


ブロック要因が有意な場合の多重比較
  投稿者:Hiroshi 2017/01/26(Thu) 17:22 No. 22260
青木先生、皆様
いつもお世話になっています。

今回はブロック間差が有意になった場合の多重比較について教えてください。データや計算過程は添付ファイルの通りです。

作物品種A〜Iについて生育量の違いがあるか多重比較しようとしています。地力差を考慮してブロックを設定し、ファイル左側にある生データについて繰り返しのない2元配置分散分析を行った結果、ブロック間差が有意でした(ファイル中央下側にある分散分析表、ブロックはP=0.0375)。そこでブロックの影響を除去した上で多重比較をすべきだと考え、ファイルの黄色に示した数値を用いて各データを補正し、ファイル右側にあるような補正済みデータセットを作成しました。例えば、品種AのBlock1の生データから−1.9(全体平均とBlock1の平均との差)を引いて50.3に補正しました。これを使って多重比較を行えば良いのでしょうか?

ご教示いただければありがたいです。
よろしくお願いいたします。


Re: ブロック要因が有意な場合の多重比較
  投稿者:太郎 2017/01/27(Fri) 09:24 No. 22261
 そのような補正は通常は行いません。
分散分析表にある誤差分散を用いて、各品種間の平均値の対比較(多重比較)を行うのが一般的です。
 古い本ですが、実験計画法(奥野・芳賀)が参考になります。
 

Re: ブロック要因が有意な場合の多重比較
  投稿者:Hiroshi 2017/01/27(Fri) 13:08 No. 22262
太郎様

早速のご教示、ありがとうございました。
ご紹介の本に当たってみます。


R findInterval()関数
  投稿者:明石 2017/01/23(Mon) 18:22 No. 22256
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

Rプログラムについて、ご教示いただきたいことが出てきました。

何卒どうぞ、よろしくお願いいたします。

ーーー

年齢(integer)のベクトルデータを読み込んで、
以下のような対応表に基づいて、
年齢(integer) → 年代(factor)にデータ変換しようと思い、
この掲示板で紹介されていました便利な関数 findInterval を適用しました。

#----------------------------------------
# 年齢 ⇒ 年代(10階級)
#----------------------------------------
#  -19   "01:-19"
# 20-24   "02:20-24"
# 25-29   "03:25-29"
# 30-34   "04:30-34"
# 35-39   "05:35-39"
# 40-44   "06:40-44"
# 45-49   "07:45-49"
# 50-54   "08:50-54"
# 55-59   "09:55-59"
# 60-    "10:60-"
#----------------------------------------

作成したRコードは、以下です。

jun <- findInterval(年齢, c(0, 20, 25, 30, 35, 40, 45, 50 ,55, 60))
cut <- c("01:-19", "02:20-24", "03:25-29", "04:30-34", "05:35-39", "06:40-44", "07:45-49", "08:50-54", "09:55-59", "10:60-")
年代 <- factor(jun, label=cut)

色々なデータへ適用したところ、
例えば、年齢が「19 以下」の区間のデータがないなどのように、
findInterval で指定した区間にデータがない場合には、エラーになるようです。

私のプログラムにミスがあるのでしょうか?

青木先生にご教示いただきたいことは、以下の2点です。
何卒どうぞ、よろしくお願いいたします。

(1) 指定した区間にデータがない場合に、エラーになる事象は、
  関数の仕様として、正しい振る舞いでしょうか?

(2) 上記エラーを回避する対処方法があるでしょうか?

ご教示いただければ、大変に助かります。
何卒どうぞ、よろしくお願いいたします。

Re: R findInterval()関数
  投稿者:青木繁伸 2017/01/24(Tue) 10:10 No. 22258
findInterval でエラーになっているのではありません。
set.seed(12345)
年齢 <- sample(18:63, 10, replace=TRUE)
年齢
jun <- findInterval(年齢, c(0, 20, 25, 30, 35, 40, 45, 50 ,55, 60))
cut <- c("01:-19", "02:20-24", "03:25-29", "04:30-34", "05:35-39", "06:40-44", "07:45-49", "08:50-54", "09:55-59", "10:60-")
年代 <- factor(jun, label=cut)
では,
> jun
[1] 8 9 8 9 5 3 4 6 8 10
> table(jun)
jun
3 4 5 6 8 9 10
1 1 1 1 3 2 1
となり,factor は 7 種類の値しか取らないのに cut が 10 種のラベルであるため,「 factor(jun, label = cut) でエラー: 無効な 'labels' です; 長さ 10 は 1 または 7 であるべきです」のエラーになるのです。
対処法は,factr の実際の値に対してlabels を用意するということになるでしょう。
#set.seed(12345)
年齢 <- sample(18:63, 10, replace=TRUE)
年齢
int = 5 # 階級幅
from <- floor(min(年齢)/int)*int
to <- ceiling(max(年齢)/int)*int
vec <- seq(from, to, by=int)
vec
jun <- findInterval(年齢, vec)
index <- as.integer(names(table(jun))) # 階級が連続しているとは限らない
labels <- NULL
for (i in seq_along(table(jun))) {
value <- (index[i]-1)*int+from
labels <- c(labels, sprintf("%02d:%d-%d", i, value, value+int-1))
}
labels
年代 <- factor(jun, labels=labels)
年代
data.frame(年齢, jun, 年代)

【御礼】 Re: R findInterval()関数
  投稿者:明石 2017/01/24(Tue) 12:57 No. 22259
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

苦慮していましたので、有難いご教示をいただき、大変に助かりました。

「factor の実際の値に対してlabels を用意する」ことで対処する方向で検討します。

いつもありがとうございます。


正規分布以外のデータの母分散の区間推定法は?
  投稿者:こびぃ 2017/01/19(Thu) 18:16 No. 22253
母平均を区間推定する場合,生データが正規分布に従わなくとも,中心極限定理のおかげで,
ある程度のサンプルサイズ(n=10個〜20個以上)があれば母平均の区間推定ができるといわれます。

実際,正規分布以外の乱数を使ったシミュレーションにより,そのことを確認しました。

一方,母分散を区間推定したい場合は,
生データが正規分布に従うことがわかっている場合,カイ2乗分布を利用します。

しかし,生データがどのような分布に従うかわからない(少なくとも正規分布ではなさそうな)場合に,
生データの母分散を区間推定する方法というのは世の中にあるのでしょうか?

正規分布よりも両裾の厚い分布(ジョンソンSU分布)を想定して母分散の推定をシミュレーションしてみたところ,
nを増やせば増やすほど的中率が低下する(信頼区間内に母分散を含まないケースが増える)という結果になりました。

これは,カイ2乗分布に従わないことに由来するのでしょう。

自分なりに調べてみたのですが,
正規分布を想定しないケースについての母分散の区間推定のやりかたが見つけられませんでしたので,
こちらの掲示板に質問させていただく次第です。

以上,ご教示いただきますようお願いいたします。

Re: 正規分布以外のデータの母分散の区間推定法は?
  投稿者:鈴木康弘 2017/01/23(Mon) 08:01 No. 22254
 カイ2乗分布が正規分布に従う標本の2乗和の分布であるように、そのジョンソンSU分布というものの2乗和の分布に従うのでは、と思うのですが..全然自信ありません。

Re: 正規分布以外のデータの母分散の区間推定法は?
  投稿者:青木繁伸 2017/01/23(Mon) 11:13 No. 22255
Efron , Stein : The Jackknife Estimate of Variance - Project Euclid
https://projecteuclid.org/euclid.aos/1176345462

Bradley Efron (1979). “Bootstrap Methods: Another Look at the Jackknife”. The Annals of Statistics 7 (1): 1–26.

Efron, B. (1981). Nonparametric estimates of standard error: The jackknife, the bootstrap and other methods. Biometrika, 68, 589-599.

Efron, B. (1982). The jackknife, the bootstrap, and other resampling plans. Society of Industrial and Applied Mathematics CBMS-NSF Monographs, 38.

簡単な関数を書いてみた
f = function(x, sig=0.05, trial=1000) {
z = double(trial)
for (i in seq_len(trial)) {
z[i] = var(sample(x, replace=TRUE))
}
quantile(z, c(sig/2, 1-sig/2))
}
これを使って,一様乱数の分散(1/12)を求めてみる
trial = 1000
result = double(trial)
for (i in seq_len(trial)) {
x = runif(50)
lr = f(x, trial=200)
result[i] = lr[1] <= 1/12 && 1/12 <= lr[2]
}
mean(result)
0.91〜0.94 と若干狭いようだけど。

最近では,ベイズ推定でもできるのかも知れない(私は知らない)

Re: 正規分布以外のデータの母分散の区間推定法は?
  投稿者:こびぃ 2017/01/23(Mon) 19:55 No. 22257
鈴木様 ご返信ありがとうございます。

ご指摘のとおり,(正規分布でなくとも)分布が具体的にわかっていれば,その分布の2乗和の分布を調べれば解決します。

しかし,現実には分布が具体的にはわからないことが多いので,その場合どうすればよいかと思った次第です。

青木先生 ご返信ありがとうございます。

具体的な文献やプログラムのご提示ありがとうございます。

実は,私もBootstrap というキーワードにこの2〜3日でようやくたどり着いたところでした。

ご紹介いただいたものを含め,Bootstrap法とかJackknife法というものを調べてみようと思います。

鈴木様,青木先生,今後ともよろしくお願いいたします。


AICの解釈
  投稿者:彷徨う初心者 2017/01/09(Mon) 21:58 No. 22248
いつも拝見させていただいています。
質問があり、書き込みさせていただきました。

A.B.Cの3つのモデルの生存曲線の比較を行っています。(いずれのモデルもスコアは1.2.3の3値をとります。)

A、Bは各群で有意に層別化されましたが、Cは1 vs 2において有意な差は認められませんでした(1 vs 3, 2vs 3は有意差あり)。

それにもかかわらず、AICはCが最も低値となります。
Cは統計学的には優れたモデルと解釈すべきなのでしょうか?

視覚的には層別化がされていないので、当然ながら良いモデルには見えません。
なぜこういった現象が起こるのでしょうか?

どなたかご教授いただければ幸いです。

p.s) この場合評価方法はC-index等にすべきでしょうか?

Re: AICの解釈
  投稿者:鈴木康弘 2017/01/12(Thu) 07:07 No. 22249
念のためですが、各モデルの説明変数の数は同じなんですよね?

Re: AICの解釈
  投稿者:彷徨う初心者 2017/01/15(Sun) 21:42 No. 22250
ありがとうございます。
サーバーの不調かアクセス出来ずスマホで返事いたします。

説明変数の数は同じです。

よろしくお願い申し上げます。

Re: AICの解釈
  投稿者:鈴木康弘 2017/01/17(Tue) 06:53 No. 22251
 AICが低値なモデルとは、予測に使えるモデルということで、群の差を出す目的とは一致しないということなんでしょうね。

 どのモデルを使うかは、目的に応じてお好きなように、ということになるのでしょう。

Re: AICの解釈
  投稿者:彷徨う初心者 2017/01/17(Tue) 11:24 No. 22252
ご返事ありがとうございます。
やはり、皆様もそういう理解ですね。

生存曲線はやはり、層別化がされている方がよいと思いますので、統計の限界かもしれません。


交絡因子の取り扱い
  投稿者:小林 2016/12/21(Wed) 10:14 No. 22243
青木先生

度々のご質問で申し訳ございません。
小林と申します。重回帰分析の質問に対して、ご丁寧にご説明して頂きまして、御礼申し上げます。
交絡因子の取り扱いについてご教授をお願い申し上げます。

交絡因子となる要因をサンプルから除外して解析する場合と、サンプルから除外しないで多変量解析時に説明因子に組み込んで解析する場合とでは結果は同じになるという認識で良いのでしょうか?

具体的には、脳卒中患者に投与している薬剤Aが脳卒中患者の認知機能を低下させているかも知れないという仮説を立証する場合に、
母数(脳卒中患者)から認知症を合併している患者を除外する場合と
認知症合併患者を除外せず、多変量解析時に認知症合併患者を説明因子として組み込む場合とでは、結果は同じになりますでしょうか?

基本的なご質問で誠に恐縮でございますが、ご教授頂ければ幸いです。
何卒、宜しくお願い申し上げます。

Re: 交絡因子の取り扱い
  投稿者:青木繁伸 2016/12/21(Wed) 12:26 No. 22244
同じにはならないでしょう。

なお,あなたの使った「母数」という用語は,用法を間違えています。

Re: 交絡因子の取り扱い
  投稿者:小林 2016/12/21(Wed) 15:13 No. 22245
青木先生

早々のお返事ありがとうございます。
両者の違いをどのように解釈したら良いかご教示願いますでしょうか?

宜しくお願い申し上げます。

Re: 交絡因子の取り扱い
  投稿者:青木繁伸 2016/12/21(Wed) 22:12 No. 22246
単変量解析の結果を積み上げても多変量解析にはならないというのは,当たり前の話。

データを層別して分析したのと,群変数も含めて分析した結果が同じ訳がないというのは当たり前。

「どのように解釈」というか,当たり前すぎて,話にならない。

解析対象を限定するというのは,「群を表す(群の特性を表す)変数を除外する」ということと同じ。必要な変数を除くというのと同じ意味を持つ。

多変量解析というのは,「変数セット」での話。
「変数セット」が違うだけで,結果が違うのは当たり前の話。
じゃあ,どの結果が正しいの?というなら,「必要な変数を揃えて分析した結果が正しいだろう」ということ。
更に言えば,「用意された変数が十分なのか?それらの変数の情報を十分に利用した分析なのか?」という質問に答えられないと行けないだろうというのが教科書的な答え。

教科書的な分析ができるかどうかは別のおはなし。

実験化学的に考えて見れば分かると思うのだけど,母集団に,A,B の二つの属性を持つ個体があるとき,その属性も関与する多変量事象(まあ,大ざっぱにいって疾患とか)についてデータをとって分析するとき,A, B のデータをどのような割合で取るか?ということを考えましょう。A:B は 1:1 でよいのか?母集団で A:B が 1:1 でないなら,それはダメでしょう。A:B = 1:10 とか A:B = 1000:1 なんていうデータを取ったらどんな悲惨な結果になるか。だからといって,A だけのデータを取ったら,それは A:B = 1:0 なんだから,よっぽど変な話になりそうだというのは誰にもわかるのでは?

また,「(疫学)調査には対照群が必要」というのもよく言われる話。「ある疾患に,ある要因が関係ある」という話をするときに,「ある疾患の患者」だけを調査して,「ある要因」が関係あるかないかは,いえないですね。

Re: 交絡因子の取り扱い
  投稿者:小林 2016/12/22(Thu) 11:45 No. 22247
青木先生

詳細な回答ありがとうございます。
申し訳ありません。私の言葉足らずでして、
どちらのやり方が正当なのかをお伺いしたかったまでです。
理解できました。
ありがとうございました。


重回帰分析の説明変数について
  投稿者:小林 2016/12/19(Mon) 10:48 No. 22240
青木先生

はじめまして。小林と申します。
こちらの掲示版で勉強させて頂いております。
今回、基本的な質問で恐縮ですが、重回帰分析の説明変数についてご質問させて頂きたく存じます。

重回帰分析の説明変数には、名義変数をダミー変数化すれば回帰式に組み込めるという認識でいたのですが、一方では、名義変数を回帰式に組み込む事はできないと仰っている先生もおられます。その理由としては、ダミー変数は、一般的には「0、1データ」として取り扱うと思いますが、これは、「0、100」として取り扱っても良い訳で、この値によって目的変数Yの値が変動するためとの事とです。

この先生の仰っている事はわかるのですが、どの統計の参考書を読んでも名義変数を回帰式に組み込む事はできないと記載されているものはありません。

重回帰分析を行うにあたり名義変数の取扱いについてどのように解釈したら良いのかご教示して頂けると幸いです。

何卒宜しくお願い申し上げます。

Re: 重回帰分析の説明変数について
  投稿者:青木繁伸 2016/12/19(Mon) 12:46 No. 22241
0/1 でなく 0/100 にしたら,そのダミー変数に対する変化いい係数が 1/100 になるだけです。目的変数 Y の値は変わりません。

以下では R でやってみますが,あなたの使える統計解析プログラムでもやってみてください
> d1 = c(0,1,1,0,0,1,1)
> d2 = c(0,0,0,1,1,0,0)
> y = c(2,5,3,6,8,9,10)
> data.frame(d1, d2, y)
試しに,このようなデータ(ダミー変数が2個)で重回帰分析をしてみると
  d1 d2  y
1 0 0 2
2 1 0 5
3 1 0 3
4 0 1 6
5 0 1 8
6 1 0 9
7 1 0 10

> a1 = lm(y~d1+d2)
> summary(a1)
重回帰分析の結果
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 2.000 2.947 0.679 0.535
d1 4.750 3.295 1.441 0.223
d2 5.000 3.610 1.385 0.238

Residual standard error: 2.947 on 4 degrees of freedom
Multiple R-squared: 0.3665, Adjusted R-squared: 0.0498
F-statistic: 1.157 on 2 and 4 DF, p-value: 0.4013
予測値は,以下の通り
> predict(a1)
1 2 3 4 5 6 7
2.00 6.75 6.75 7.00 7.00 6.75 6.75


ダミー変数を 0/1 ではなく 0/100 にしてみる
目的変数は当然同じ
> e1 = d1*100
> e2 = d2*100
> data.frame(e1, e2, y)
つまりこんなデータを重回帰分析
   e1  e2  y
1 0 0 2
2 100 0 5
3 100 0 3
4 0 100 6
5 0 100 8
6 100 0 9
7 100 0 10

> a2 = lm(y~e1+e2)
> summary(a2)
分析結果は以下の通り
偏回帰係数(Estimate)は前の結果の 1/100 になっている
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 2.00000 2.94746 0.679 0.535
e1 0.04750 0.03295 1.441 0.223
e2 0.05000 0.03610 1.385 0.238

Residual standard error: 2.947 on 4 degrees of freedom
Multiple R-squared: 0.3665, Adjusted R-squared: 0.0498
F-statistic: 1.157 on 2 and 4 DF, p-value: 0.4013
予測値は以下の通り。
前の重回帰分析での予測値と同じ(同じになるように偏回帰係数を調整したことになる)
> predict(a2)
1 2 3 4 5 6 7
2.00 6.75 6.75 7.00 7.00 6.75 6.75

ここでは話をわかりやすくするために 0/1 を 0/100 にしたが,二値変数であればなんでもよい。0/3.14 でも 1:110 でも 3.14/123.45 でも予測値には何の違いもない。
0/1 にしているのは,その要因があるかないかということに対応するので結果の解釈が容易になるというだけの話。

Re: 重回帰分析の説明変数について
  投稿者:小林 2016/12/21(Wed) 10:02 No. 22242
青木先生

早速のお返事ありがとうございます。
自分の統計ソフトでも確認致しました。
目的変数は変わりはなく、偏回帰係数のみが変わるのですね。
理解できました。

ご対応頂きありがとうございました。


因子分析の信頼性と妥当性について
  投稿者:さくら 2016/12/07(Wed) 19:01 No. 22237
先生の掲示板でいつも色々と学ばせて頂いております。

これまでは因子分析で抽出された因子の信頼性評価をクロンバックαに
より求めてきました。

先行研究に、合成信頼性を使った信頼性評価と、平均分散抽出に
よる妥当性が必要と記載の上、検証されています。現在、アンケートに
より回収したデータを用いて因子分析を行っている最中です。

平均分散抽出とはどのような意味で、どのように値を求めたらよいので
しょうか。また合成信頼性についても意味するところと値の求め方が
分からずに困っております。

今後、このあたりについてどのような文献を用いて学んでいったら
良いのか、目星がつけられません。

ご教示頂ければ幸いです。よろしくお願い申し上げます。

Re: 因子分析の信頼性と妥当性について
  投稿者:鈴木康弘 2016/12/09(Fri) 17:36 No. 22238
 日本語のサイトではこれというのが見つかりませんね。
"average variance extracted" と "composite reliability" で探すと

https://www.researchgate.net/post/Can_anyone_tell_me_how_to_calculate_average_variance_extracted_AVE_and_composite_reliability_CR_of_a_single_latent_variable_with_5_indicators

 というのが当たります。英語ですけど。

Re: 因子分析の信頼性と妥当性について
  投稿者:さくら 2016/12/15(Thu) 15:23 No. 22239
鈴木様

早々にサイトをご紹介いただきありがとうございました。
キーワードまで示していただき、今後の勉強の参考になりました。

ご紹介いただいたサイトで、まずは勉強してみます。

改めてお礼申し上げます。


R カラム数・カラム名が異なる複数ファイルの、行方向での結合
  投稿者:明石 2016/12/01(Thu) 22:10 No. 22227
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

Rプログラムについて、ご教示いただきたいことが出てきました。

何卒どうぞ、よろしくお願いいたします。

ーーー

データフレームの行方向での結合 rbind()関数を適用できる前提条件は、
カラム数・カラム名の並びが同一であることは理解しております。

カラム数・カラム名が異なる複数ファイルの、行方向での結合について
ご教示をいただければ大変に助かります。

結合対象のファイル名は、fllenames[n]のベクトルに格納されているとします。

素人考えですが、以下のようにやればよいのかと、思いを巡らせている段階です。

・まず、fllenames[1]のファイル名で読み込み、データフレームdat1に格納します。

・次に、fllenames[2]のファイル名で読み込み、データフレームdat2に格納します。

・dat1のカラム名と、dat2のカラム名とを比較して、
 rbindができるように、dat1とdat2のカラム名を整形してあわせる(※)。

・これを、逐次的に、fllenames[n]まで、すべてのファイル名について行う

言葉ではさらりと書きましたが、※ の箇所の具体的な手順を書き出すことができません。

もっと良い方法があるように思います。

もし可能ならば、具体的なRプログラムをお示ししていただき、
ご教示をいただければ大変に助かります。

何卒どうぞ、よろしくお願い申し上げます。

失礼いたします。

Re: R カラム数・カラム名が異なる複数ファイルの、行方向での結合
  投稿者:青木繁伸 2016/12/02(Fri) 09:05 No. 22228
具体的によくわからないけど,merge でできそうに思いますが?

Re: R カラム数・カラム名が異なる複数ファイルの、行方向での結合
  投稿者:明石 2016/12/02(Fri) 09:56 No. 22229
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

早々にご教示をいただきまして、誠にありがとうございました。

「転置させて、カラム名をキーにしてmergeをする」という意味と理解しました。

これだけでも貴重な手がかりです。
これを手がかりにして、検討をしてみます。

分かりにくい、つまらない質問で、大変に申し訳ございませんでした。
以降、気をつけます。
失礼いたします。

Re: R カラム数・カラム名が異なる複数ファイルの、行方向での結合
  投稿者:青木繁伸 2016/12/02(Fri) 11:15 No. 22230
具体例を示してくれればよいのだけど?
merge は(大部分が)同じ個体についてのデータフレームをマージするのが目的。
「カラム数・カラム名が異なる複数ファイル」というのは,全く別の個体のデータですか(普通はそんなのは考えられないけど)。
もしカラム名の重複がないデータなら,両者に同じカラム名(以下では id)で値の重複もないカラムを追加して,all=TRUE で merge すればよい。
a = data.frame(id=1:5, x=1:5*3, y=(1:5)^2)
a
b = data.frame(id=11:20, z=runif(10))
b
merge(a, b, all=TRUE)
実行例
> a = data.frame(id=1:5, x=1:5*3, y=(1:5)^2)
> a
id x y
1 1 3 1
2 2 6 4
3 3 9 9
4 4 12 16
5 5 15 25
> b = data.frame(id=11:20, z=runif(10))
> b
id z
1 11 0.45929778
2 12 0.03272818
3 13 0.51413244
4 14 0.07282785
5 15 0.78663924
6 16 0.12525008
7 17 0.87019710
8 18 0.37075708
9 19 0.13116714
10 20 0.44714621
> merge(a, b, all=TRUE)
id x y z
1 1 3 1 NA
2 2 6 4 NA
3 3 9 9 NA
4 4 12 16 NA
5 5 15 25 NA
6 11 NA NA 0.45929778
7 12 NA NA 0.03272818
8 13 NA NA 0.51413244
9 14 NA NA 0.07282785
10 15 NA NA 0.78663924
11 16 NA NA 0.12525008
12 17 NA NA 0.87019710
13 18 NA NA 0.37075708
14 19 NA NA 0.13116714
15 20 NA NA 0.44714621


Re: R カラム数・カラム名が異なる複数ファイルの、行方向での結合
  投稿者:明石 2016/12/03(Sat) 09:45 No. 22231
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

ご親身にご相談にのってくださり、誠にありがとうございます。
改めて御礼を申し上げます。

やりたいことをどのようにお伝えすればよいのか、色々と考えており、
ご返信が遅くなり、申し訳ございませんでした。

青木先生がお書きになられたことは、よく理解できました。
ありがとうございます。

私がやりたいと思うことは、
個体〜観測変数 というデータの結合ではなく、
データベースから出力した帳票類のデータファイルの結合です。

入力ファイルが3つからなる簡単な例示を、添付画像でお示しをします。
ファイル1の項目は、v1、v2、v3
ファイル2の項目は、v1、v2、v4
ファイル3の項目は、v2、v3、v5

これら3つのファイルを行方向に結合したいのですが、
カラム名が一致しないために、単純にはrbindが適用できません。これが悩みです。

共通する項目については、それを利用し、
共通しない項目については、右側の列に追加しながら、
行方向に結合します。

データがない箇所は、NAで埋めます。

結果として、
行方向の数は、入力ファイルの行数の和
列方向の数は、入力ファイルのカラム名の和集合
となります。

これがやりたいことです。

先の投稿では、説明不足で申し訳ございませんでした。
どうぞ、よろしくお願い申し上げます。


Re: R カラム数・カラム名が異なる複数ファイルの、行方向での結合
  投稿者:青木繁伸 2016/12/03(Sat) 13:18 No. 22232
なので,merge を使えばできるでしょう?
> file1 = data.frame(v1=rep(1,5), v2=rep(1,5), v3=rep(1,5))
> file2 = data.frame(v1=rep(2,5), v2=rep(2,5), v4=rep(2,5))
> file3 = data.frame(v2=rep(3,5), v3=rep(3,5), v5=rep(3,5))

> file1$id = 1:5 # ユニークな値を取る変数を付加する
> file2$id = 1:5 + 10 # 不要な場合もあるが
> file3$id = 1:5 + 20 # 間違いを防ぐためには必須

> file1
v1 v2 v3 id
1 1 1 1 1
2 1 1 1 2
3 1 1 1 3
4 1 1 1 4
5 1 1 1 5
> file2
v1 v2 v4 id
1 2 2 2 11
2 2 2 2 12
3 2 2 2 13
4 2 2 2 14
5 2 2 2 15
> file3
v2 v3 v5 id
1 3 3 3 21
2 3 3 3 22
3 3 3 3 23
4 3 3 3 24
5 3 3 3 25

> a = merge(merge(file1, file2, all=TRUE), file3, all=TRUE)
> a = a[,c("v1", "v2", "v3", "v4", "v5")] # 列順を決める必要があるならこれを追加
> a
v1 v2 v3 v4 v5
1 1 1 1 NA NA
2 1 1 1 NA NA
3 1 1 1 NA NA
4 1 1 1 NA NA
5 1 1 1 NA NA
6 2 2 NA 2 NA
7 2 2 NA 2 NA
8 2 2 NA 2 NA
9 2 2 NA 2 NA
10 2 2 NA 2 NA
11 NA 3 3 NA 3
12 NA 3 3 NA 3
13 NA 3 3 NA 3
14 NA 3 3 NA 3
15 NA 3 3 NA 3
複数のデータフレームのマージは,
a = read.table(filenames[1])
for (i in 2:n) {
a = merge(a, read.table(filenames[i]), all=TRUE)
}
のようにやる。

【御礼】 Re: R カラム数・カラム名が異なる複数ファイルの、行方向での結合
  投稿者:明石 2016/12/03(Sat) 17:02 No. 22233
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

青木先生が「mergeでできる」とおっしゃった意味を、納得し、理解できました。

今回も助けていただきました。
心から御礼を申し上げます。
ありがとうほざいました。

失礼いたします。

Re: R カラム数・カラム名が異なる複数ファイルの、行方向での結合
  投稿者:明石 2016/12/04(Sun) 08:44 No. 22234
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

おさらいをしていて、追加質問です。
どうぞ、よろしくお願いいたします。

私は、データベース操作言語のSQLを利用していますので、
Rのmergeも、その理解でおりました。

結合の種類(内部、左外部、右外部、完全外部)に応じた指定方法がありますが、
結合キーは必須だと思っていました。

青木先生からご教示いただきました方法では、結合キーの指定がないことに、
驚きました。

Rのmergeは、SQLの結合だと思っていた私には、
このような使い方があることを、まったく知りませんでした。

追加質問です。

> file1$id = 1:5 # ユニークな値を取る変数を付加する
> file2$id = 1:5 + 10 # 不要な場合もあるが
> file3$id = 1:5 + 20 # 間違いを防ぐためには必須

idを付加する場合、しない場合について、
私なりに色々と実験をしましたが、この意味合いが理解できません。

ご説明をいただけますと、大変に助かります。
どうぞ、よろしくお願いいたします。

色々とありがとうございます。

Re: R カラム数・カラム名が異なる複数ファイルの、行方向での結合
  投稿者:青木繁伸 2016/12/04(Sun) 11:18 No. 22235
結合キーは,R が自動で指定します
両者に共通の列です
今のような場合には複数のファイルを指定することもあり,id 列を付加することで,面倒な処理を避けることができるでしょう
id 列を付加しないと,思いがけない結果になりびっくりすることになります
以下の例を吟味するとよいでしょう(file2 が前の例と異なります)
file1 = data.frame(v1=rep(1,5), v2=rep(1,5), v3=rep(1,5))
# file2 = data.frame(v1=rep(2,5), v2=rep(2,5), v4=rep(2,5))
file2 = data.frame(v1=rep(1,5), v2=rep(1,5), v4=rep(2,5))
file3 = data.frame(v2=rep(3,5), v3=rep(3,5), v5=rep(3,5))
file1
file2
file3
a = merge(merge(file1, file2, all=TRUE), file3, all=TRUE)
a = a[,c("v1", "v2", "v3", "v4", "v5")] # 列順を決める必要があるならこれを追加
a

Special Thanks(R カラム数・カラム名が異なる複数ファイルの、行方向での結合)
  投稿者:明石 2016/12/04(Sun) 12:13 No. 22236
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

何度も何度もご丁寧にご教示をいただきまして、誠にありがとうございました。

おかげさまで理解できました。

新しい追加拡張パッケージが沢山ありますが、
baseの基本機能を正しく理解することで、対応できることが多いなと、
改めて、Rはよく設計されていると感心しました。

心から御礼を申し上げます。
失礼いたします。


カイ2乗検定について
  投稿者:心理系院生 2016/11/30(Wed) 12:56 No. 22224
2群の割合に違いがあるかを検討するために、カイ2乗検定を行うのがよいと考えました。
しかし、2群の標本数が大きく異なっており、このような場合にカイ2乗検定を行ってよいものかわかりませんでした。
一方は10000対7000ほどで、もう一方は60対10程度です。
このような違いがあってもカイ2乗検定を行ってもよいのでしょうか。
もし他に適切な手法があれば、ご教授いただけましたら幸いです。
お忙しい所恐れ入りますが、何卒よろしくお願いいたします。

Re: カイ2乗検定について
  投稿者:かい 2016/12/01(Thu) 13:14 No. 22225
サンプルサイズが異なっても検定に問題ありませんが,10000対7000の方は十分精度が高い(区間推定の範囲が非常に狭い)ので,実質的には10000対7000の方の点推定値と60対10の方の比率を比較しているのとあまり変わらない結果になります.

Re: カイ2乗検定について
  投稿者:心理系院生 2016/12/01(Thu) 14:54 No. 22226
かい様

お返事ありがとうございます。
よく理解できました。
問題はないということで安心しました。
このまま研究を進めていきたいと思います。


尺度開発
  投稿者:kashio 2016/11/15(Tue) 14:19 No. 22211
コミュニケーションに関する尺度を開発中です。
いくつかの書籍を参考に因子分析をしました。累積寄与率を示したかったし、スクリュープロットで4因子にかくていし、下位尺度や信頼性の結果も納得できたので投稿したところ査読者から以下のコメントをもらいました。今時 が引っかかり、納得も行かず、相関関係が認められない4因子に分けられる確証のところの意味もわからないのです。
素直に、解析をし直すべきでしょうか?

主因子法とバリマックス回転を用いた因子分析をされていますが、理由が見当たりません。言えなかった状況9カテゴリーと、言わなければよかった状況4カテゴリーから作成された尺度が、必ず相関関係が認められない4因子に分けられるという確証が無い限りは、今時は最尤法とプロマックス回転などの斜行回転を用いるべきであると考えます。

Re: 尺度開発
  投稿者:ポン太 2016/11/15(Tue) 15:05 No. 22212
どのような書籍を参考にされたのかわかりませんが、査読者の指摘の通りだと思います。

主因子法は別にして、実質科学からして、相関関係をゼロと仮定するパリマクスが理由もなく
使われていたら、因子分析を理解されている査読者なら、??となりますよ。

納得いかないのであれば、編集委員会に伝えたらいいでしょう。しかし今回はおそらく却下されてしまいますよ。

良書を選んで読まれてはいかがでしょうか。

Re: 尺度開発
  投稿者:kashio 2016/11/15(Tue) 15:28 No. 22213
ありがとうございます。良書をご紹介いただければと思います。
プロマックス回転に挑戦します。この査読者は、これだけですが、
もう一人の査読者は以下のように言っています。共分散は次の論文でと思っていたのですが、プロマックスに挑戦してもダメでしょうか?

α係数も、0.6台が3因子と、正直高いとは言えません。そのことも考えると、確証的因子分析で本尺度項目での因子構造の確認を行うことまでしたほうがよいのではないでしょうか。

Re: 尺度開発
  投稿者:ポン太 2016/11/15(Tue) 15:49 No. 22214
はい、確認的因子分析まで持ち込む必要があると思います。

「共分散はつぎの、、、、 ダメでしょうか」の意味が理解できません。分析には、分散共分散行列、もしくは相関行列が使われているはずなのですが。

いずれにせよ、あなたがやるべきことは反論ではなく、2人の査読者の指摘されたことだと思います。

Re: 尺度開発
  投稿者:ポン太 2016/11/15(Tue) 15:51 No. 22215
書き忘れましたが、因子間相関がゼロと言えるのであれば、反論の余地はあるでしょう。

Re: 尺度開発
  投稿者:kashio 2016/11/15(Tue) 17:27 No. 22216
ありがとうございました。確証的因子分析してあります。今回は質的研究で質問項目を起こしました。質的研究で分けたカテゴリーから項目を作っていった過程の記述を丁寧に書いてあるので紙面が足りず探索的因子分析までにしようと思いました。頑張ってみます。

お願いです。良書を紹介してください。

Re: 尺度開発
  投稿者:kashio 2016/11/15(Tue) 22:24 No. 22217
ポン太さま

やっぱり、わかったふりはやめます。今時 に引っかかっていたのは、何年か前に還暦をすぎていて、主因子法バリマックス回転で何本か論文を書いているから!でも、今時 だからではなく、因子間相関がゼロでなければ、バリマックス回転はしない。ということがわかりました。バリマックスで一番納得がいく下位尺度ができたのですがー
残念ですが、最尤法プロマックスでやってみます。
確証的因子分析はやってあるし、いわゆる当てはまりもいいようです。
でも、なぜ、確証的因子分析をしなくてはならないのかわかっていません。
反論をしたいのではなく、理解したいのです。査読者は、アルファー係数が0.6と低いので、確証的因子分析を進めていますが、アルファ係数は信頼係数だと理解しています。ご助言ください。

Re: 尺度開発
  投稿者:ポン太 2016/11/16(Wed) 08:43 No. 22218
確認的因子分析をされているので。その意味は理解されていると思ったのですが…。因子構造の単純化です。探索的および確認的因子構造を作図してみれば一目瞭然かと。作図しなくても理解はできますが。

当てはまりはいいかもしれませんが、因子間相関をゼロと仮定するか否かのところで、近年の統計学の理論が反映されていないのが問題です。ただわたしは、純粋な統計学者ではありませんので、数式などを使っての説明はできません。

アルファ係数は高いに越したことはないでしょうが、大学の成績の秀優良可くらいに考えておけばいいのでは?0.6なので可ですか。まあ、この基準を論文に書けるかは別問題になりますが。そもそもこうあるべきというルールはないと思います。アルファ係数は内的信頼性係数のことでいいと思います。青木先生のRの関数を使って、係数を高めることは可能だと思います。もちろん数値だけをみてバンバン削除するのは違いますが。

豊田秀樹氏、狩野裕氏の著書が参考になります。もちろん他にもあります。

Re: 尺度開発
  投稿者:kashio 2016/11/16(Wed) 10:38 No. 22219
ポン太様

ありがとうございます。「理解した」とは自分の言葉で語れること,と思っています。
まだ語れません。きちんと理解できるように、まず、ご紹介いただきました書籍を読みます。

Re: 尺度開発
  投稿者:青木繁伸 2016/11/16(Wed) 23:29 No. 22220
因子間に相関があることを仮定しているのがプロマックス回転であるからして,バリマックス回転はそれに含まれる(相関が0である特殊な場合)。
よって,プロマックス回転をして,ある因子間の相関はほとんど0であるといえばよいだけの話。すべての因子間相関が0に近いならば,「バリマックス回転でよいのだ」といえるわけだけど,因子間相関が0であることをことさらい言う意味がどれほどあるかということ。

まあ,昔は(コンピュータが使えなかったり(そんなばかな),斜交回転までは検討できなかったりで)主因子解・バリマックス回転を行うのが主流でもあったのだけど,現在に至っては「いまさら」なんでしょうね。査読者の「今の時点じゃ,そんなの当たり前でしょ」という,上から目線がや〜〜な感じなのは仕方ない(あなた方も,少し前までは主因子解でバリマックス回転で論文書いてたでしょ...と今更いっても仕方ない)。

確証的因子分析の必要性も,初期の頃には重視されていたとは言いがたい。

尺度が使用に耐えるためにはアルファ係数は 0.8 以上というのは教科書的ではあります。0.7 以上ならまあいいのではないかとかいう人もいたりするわけではありますが,0.6 で十分というのは,私は見聞きしたことはないですね。

統計学を用いた論文の論文たる水準も時代によって変化する(上昇する)ので,そのときどきの水準を満たすことは学術論文として必須のことではあるのでしょう。

いずれにしろ,学会誌に論文を掲載してもらうにはまずは査読者の判定を経なければならないわけで,正当な反論がない限り泣く子と地頭には勝てないので,査読者の指示に従うしかないということになるでしょう。

査読者と喧嘩をしてみるのも面白いかも知れませんね。

Re: 尺度開発
  投稿者:kashio 2016/11/18(Fri) 16:12 No. 22221
ポンタ様、青木先生ありがとうございました。
最尤法、プロマックス回転をしました。
主因子法、バリマックス回転と同じ因子、項目数でした。因子固有値が高くなりかえって良い結果でした。下位尺度ごとに相関はありました。
もう一つ査読者が 項目分析の結果を書けと言ってきました。
これが良くわかりません。分布の偏りや、標準偏差を確かめたり
しました。下位尺度ごとにクロンバックを上げるように項目を出し入れしました。
その都度テスト再テストや基準関連妥当性も見ました。

クロンバックですが、査読者は0.6しかないといっていますが、
第1因子0.68、第2因子0.75、第3因子0.68、第4因子0.68 全体0.81
でした許されませんか?隣の部屋の先生に勧められて、KMOとBartlettなるものもやって、標本妥当性なるものもみました。隣の先生は有名な先生ですけど項目分析はしらないといってます。項目分析を調べたところクロンバックがなんとかなったのでよいのではないかと思うのですが、この項目分析も今時のものですか?、

Re: 尺度開発
  投稿者:青木繁伸 2016/11/18(Fri) 21:04 No. 22222
クロンバックのアルファですが,教科書的には 0.8 というのが閾値ですが,「実用では 0.7 でもよい」とか書いてあるサイトも見かけましたが,その根拠については確認していません。

わたしも,項目分析についてはよくわかりません。

Re: 尺度開発
  投稿者:kashio 2016/11/19(Sat) 21:15 No. 22223
ありがとうございます。参考書を読んでもなぜかがわからず苦しむことがおおいです。もう少しまなんで、探索的因子分析だけでなく、なぜ確証的因子分析をしなくてはならないか理解したいです。


二つの級内相関係数の差の検定
  投稿者:NIIGATA 2016/11/09(Wed) 22:59 No. 22210
お世話になります。
唐突で大変恐縮ですが、2つの級内相関係数の差の検定なるものは存在するのでしょうか?

そして、存在するならば、それは相関係数の差の検定とはどのように異なるのでしょうか?

可能でしたらRのスクリプトも含めご教示いただけますと幸いです。


工程能力指数Cp
  投稿者:かずみ 2016/11/04(Fri) 21:32 No. 22207
「工程能力指数Cpは、規格の中心と標本平均値Xberがほぼ同じ時に求める。」となっていますが、このほぼ同じ位とはどれ位までなのですか?
Cpkは、こうしたことを考えなくて良いので、実務ではCpkを用いていますが、
Cpの場合は、この辺りが曖昧でわかりません。検定で同じとは言えないので、分散分析などのプーリングの考え方などから考えてαを25%と考えてこれを両側として、、、、すると±0.3シグマ位をほぼ同じと考えてよいのでしょうか?
このほぼ同じくらいといった事を統計的に解説した書物、論文はあるのでしょうか?

Re: 工程能力指数Cp
  投稿者:まさき君 2016/11/09(Wed) 08:19 No. 22208
品質管理上の経験則であって、理論的に説明できないとは思いますが、工程ばらつきにおいては、ロット毎に平均値は変動するが、ロット内の分散は殆ど変動しない、と考えるのが一般的で、これに基づけば、Cpにも一定の意味合いはあります。
例えば、Cpkが全然駄目な時、Cpは十分ならXberを補正する事で歩留り改善の余地があるが、Cpも駄目なら手のつけようがない、等。

Re: 工程能力指数Cp
  投稿者:かずみ 2016/11/09(Wed) 18:10 No. 22209
ありがとうございます。
やはり工程能力指数はCpkとしてXberが補正、調整可能な場合の工程能力指数がCpとかんがえるべきなのでしょうか?
であれば色々な書物、文献はそうすべきですね。


集団と個体との比較
  投稿者:さかな 2016/10/27(Thu) 16:11 No. 22194
ある個体が,集団と同じレベルといえるかどうかを調べようとしています。
既存の集団があって,その集団から1個体を抜き出してある処理をした上で,元の集団と差があるかどうかを調べるというものです。
抽象的で分かりにくいので,架空の具体例でお伝えします。

トリュフを探すための豚を80頭飼っており,その中から1頭(Aとする)を抜き出して嗅覚が良くなる薬を投与します。その後,50日間,全頭で同時に同じところでトリュフを探させる実験を毎日します。毎日異なる土地へ連れて行って,探しやすさなどの条件を変えます。日ごとに,各豚が見つけたトリュフの数を数えます。全ての豚は常に個体識別がしっかりとできる状態です。

薬の効果を確認するために,Aが有意にトリュフ探しの能力が向上したか否かを統計的に調べるにはどのようにしたら良いでしょうか。また,もし能力が向上したとしたら,どの程度向上したのか(1.4倍とか)をどのように判断すれば良いのでしょうか。

調べようとする試験区個体が1個体だけでは無理があるのかもしれませんが,この場合非常に高価な薬であるために,1頭分しか用意することができませんでした。また,Aは50回全ての日に参加していますが,その他の豚は参加しなかった日があるため,対照区の豚の頭数は日によって異なります。

Re: 集団と個体との比較
  投稿者:青木繁伸 2016/10/27(Thu) 18:17 No. 22195
実験計画に失敗しましたね
統計学的に意味のないデータです
高い薬と時間と人手が無駄になりました

Re: 集団と個体との比較
  投稿者:鈴木康弘 2016/10/28(Fri) 07:41 No. 22196
その80頭を、1番鼻のよい豚、2番目に鼻のよい豚、と並べることができるなら、薬を与えた豚が2番目までに入っていれば、5%有意水準の両側検定で有効だと、言えるかも?

Re: 集団と個体との比較
  投稿者:青木繁伸 2016/10/28(Fri) 08:45 No. 22197
豚が1頭だと,その測定値が何個得られようと,測定値が表すのは個体内変動であり,その測定値の平均値が大きくても,その豚に限ってはそうかも知れないけど,他の豚にも同じ効果があるとはいえないのではないか?に対しては反論できない。測定回数を増やしても,その豚の特性の測定精度が高くなるだけあり,効果が一般的であるかどうかについては回答できないままである。
統計学は,複数の個体について測定し,その測定値が表すのは個体間変動であり,測定値の平均が誤差変動よりも大きければ効果があるとする。効果が顕著な豚もそうでない豚もあるけど平均すれば効果がある。この場合には,サンプルサイズが大きくなれば効果の測定値の精度が高くなる(より確かな推定値が得られる)。

「薬を与えた豚が2番目までに入ってい」ても,「その豚はそうだったでしょうけど」ということに変わりはない。統計解析手法の選択の問題ではなく実験計画時点での問題。データを得たあとで分析不能があきらかになるという不都合。

Re: 集団と個体との比較
  投稿者:さかな 2016/10/31(Mon) 10:12 No. 22200
青木先生,鈴木先生

コメントいただきまして,誠にありがとうございます。
実験する豚が1頭だけではその個体に対する効果だけしかわからないので回数を重ねても同じで意味がないということ,理解いたしました。ありがとうございます。

では,この豚に対するこの薬の効果を知りたい場合には,何らかの方法で効果の有無を判断することができるでしょうか。以下のようなことを考えてみましたが,いかがなものでしょうか。

ここで,すべての豚の50回(50日)の成績(トリュフの個数)がわかっています。
例えば,50回のそれぞれにおける全体の中でのAの位置を知るために,偏差値を計算し,全ての回で偏差値66.4以上(上位5%以内)であったとしたら,有意にAの成績が良かったと考えてよろしいでしょうか。
あるいは全ての回で偏差値33.6(下位5%)〜66.4(上位5%)であった場合に,有意に成績が良かったり悪かったりすることがないといえるでしょうか。
さらに,50回のうち何回か(例えば9回)偏差値70以上のケースがあった場合には,どのように判断すればよいでしょうか。

最初と少し議論がずれてしまいますが,薬の豚に対する効果ではなく,あくまでこの個体に対する効果を見たい場合です。

Re: 集団と個体との比較
  投稿者:青木繁伸 2016/11/01(Tue) 08:24 No. 22201
その個体の一連の測定値の中の最小値が他の普通の個体のデータと一緒に外れ値の検定を行って,外れ値だといえればよいかな。

最小値を使うのは,最も控えめに言うため。平均値でもよいかも知れないけど,批判する人は批判するだろうから。

Re: 集団と個体との比較
  投稿者:さかな 2016/11/01(Tue) 09:09 No. 22202
青木先生

ありがとうございます。
少し説明不足のところがあったようです。すみません。

この実験は,毎日条件の異なる土地で実施ています。したがって,それぞれの豚が見つけたトリュフの数は1桁の時もあれば,3桁の時もあります。全体的にほかの豚が多く見つけるときはAも多く見つける傾向がみられます。したがって,各実験日のデータは独立(という用語が適切かわかりませんが)であるため,Aの測定値の最小値や平均を使うのは適切ではないのではと考えます。
つまり,各実験日のデータごとに比較をして,50日分のデータを合わせて差があるといえるのかどうかを判断したほうが良いのではないかと思います。
いかがでしょうか。

その場合には,実験日ごとに外れ値の検定をして,もちろん全実験日で外れ値またはそうでないとの結果であれば,それぞれ差がある,または無いと判断できますが,どちらの場合も含まれるときに,どのようにして判断すればよいのか,ご教授いただけたら幸甚に存じます。

Re: 集団と個体との比較
  投稿者:青木繁伸 2016/11/01(Tue) 09:43 No. 22203
> この実験は,毎日条件の異なる土地で実施ています。したがって,それぞれの豚が見つけたトリュフの数は1桁の時もあれば,3桁の時もあります。全体的にほかの豚が多く見つけるときはAも多く見つける傾向がみられます。したがって,各実験日のデータは独立(という用語が適切かわかりませんが)であるため

その一頭の豚のデータは独立ではないですね。多い日も少ない日もあるでしょうが,潜在能力は一定なのでしょうから。

しかし,まあ,ここまで「測定値」が特徴的だと,どうしようもないのではないですか。

一頭ではあるが,その一頭(他の豚も)なまじ複数回容易に測定できるから複雑になる。
複数頭からの1回のデータというのが大原則。

Re: 集団と個体との比較
  投稿者:さかな 2016/11/01(Tue) 11:46 No. 22204
すみません,例を変えてお聞きいたします。

湖にアヒルが80羽いて,そのうち1羽(Bとします)が野犬に襲われてけがをしてしまいました。嘴にけがをして,餌を食べにくい状況です。

このBを保護していますが,Bが他のアヒルと同じくらい餌を食べるのであれば,湖に返そうと思います。

そこで,1か月間毎日,Bの餌の摂取量とほかのアヒルの摂取量をそれぞれ確認しました。ただし,ここでは1カ月間での回復はないものと仮定します。

このデータによって,Bの餌の摂取量が他のアヒルと比べて有意に少ないか否かを確認することはやはり無理なのでしょうか。

Re: 集団と個体との比較
  投稿者:鈴木康弘 2016/11/03(Thu) 07:15 No. 22205
青木先生の言われるように、外れ値の検定では、いかがでしょう?

Re: 集団と個体との比較
  投稿者:さかな 2016/11/04(Fri) 13:28 No. 22206
青木先生,鈴木先生

いろいろとご教授いただきまして,誠にありがとうございました。


データ間の相関
  投稿者:初心者 2016/10/29(Sat) 13:33 No. 22198
青木先生、いつもお世話になっております。

Aという病気をもつ患者50名とBという病気をもつ患者20名に、それぞれ5つの心理検査を行いました。5つの心理検査のデータのうち、2つの心理検査間での総得点の相関や各心理検査の下位項目間での得点の相関を調べるときに用いる統計学的手法は因子分析になるのでしょうか?ちなみに結果についての仮説はありません。

初歩的な質問で大変申し訳ありませんが、どのような統計学的手法を用いたらいいのか分からず、困っています。

恐れ入りますが、どうぞご教示の程、よろしくお願い致します。

Re: データ間の相関
  投稿者:青木繁伸 2016/10/29(Sat) 18:26 No. 22199
いくつかの変数間の相関関係を見るのなら相関係数と言うことになるのでしょうけど
因子分析をするにしては,サンプルサイズが小さいし
A,B の二疾患をもつ患者をひとまとめにして分析してもよいかということもありますね
初心者とはいえ,いろいろなことを知っていないと適切な分析はできませんねえ
(そもそも,心理検査って,疾患ごとに特徴的な所見が得られたりするんでしょうかね)

[1] [2] [3] [4]

- J o y f u l  N o t e -
Modified by i s s o