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

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

[トップに戻る] [利用上の注意] [ワード検索] [過去ログ] [統計学自習ノート] [管理用]
latest article: No. 22369, 2017/05/27(Sat) 21:14
おなまえ
タイトル
コメント
ＵＲＬ
添付ファイル
暗証キー (英数字で8文字以内)
自分の投稿記事を修正・削除する場合 --> 処理 記事No 暗証キー

ベクトルが作れない
　　投稿者：課長 2017/05/26(Fri) 05:02 No. 22358
Rでベクトルを作ろうと思ったら、次のようなメッセージが出てしまいます。

> x<-c(1,1,1,1,1,1)
c(1, 1, 1, 1, 1, 1) でエラー:
6 個の引数が 'exp' に渡されましたが、1 が必要とされています

どうして急にベクトルが作れなくなったのか、ご教示いただければ幸いです。
Re: ベクトルが作れない
　　投稿者：青木繁伸 2017/05/26(Fri) 11:35 No. 22359
どこかで，c という関数を再定義してしまったのではないでしょうか。

コンソールに c と入力して
> c
function (...)  .Primitive("c")
以外のものが表示されるようならば，確認の上で，
> rm(c)
と入力して，その再定義されたものを消去すればよいでしょう。
Re: ベクトルが作れない
　　投稿者：課長 2017/05/26(Fri) 20:53 No. 22363
青木先生のおっしゃる通りでした。なぜかCに関数が定義されていました。そしてご指導通りにしたらなおりました。ありがとうございます
Re: ベクトルが作れない
　　投稿者：課長 2017/05/26(Fri) 22:49 No. 22364
青木先生の言われるとおりにしたらなおったのですが、Rを再起動するとまたｃに謎の関数が定義されてしまいます。これはどうしたら直るかご教示いただけませんか。

> c
function (x) .Primitive("exp")

が毎回出ます。
Re: ベクトルが作れない
　　投稿者：青木繁伸 2017/05/27(Sat) 07:00 No. 22365
R を起動してすぐに rm(c) をしたあと，quit(save="yes") でワークスペースを保存して R を終了してください。
そのあと R を再起動すると c の定義は永久になくなると思います。
もしそれでもなくならないなら，実行されるファイルのどこかに c <- exp というような c を再定義する箇所がないか捜してください。
なお，蛇足ですが，c を exp のような関数ではなく定数を代入しても，c 関数は置き換えられません。
> c <- 1
> c(9, 3, 1)
[1] 9 3 1
Re: ベクトルが作れない
　　投稿者：課長 2017/05/27(Sat) 21:14 No. 22369
ありがとうございます、青木先生！
おっしゃるとおりにいたしましたところ、直りました。
また、勉強にもなりました。
心より感謝いたします。
今後ともよろしくお願いいたします・


Ｒ　レイアウトが異なるテキストファイルの読み込み
　　投稿者：明石 2017/05/26(Fri) 14:33 No. 22360
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

Ｒについて、お聞きしたいことがでてきました。
ご教示をいただければ、大変に助かります。
どうぞよろしくお願いいたします。

ーーー

添付の画像ファイルをご覧ください。

読み込むテキストファイルは、
１行目
２行目〜
と列の長さが異なることから、ファイルの読み込みで苦慮しています。

データの区切りは、半角空白文字です。

1行目の２つの数字（300000 200）は、２行目以降の表のサイズを表しています。
題材によりサイズは異なります。

２行目以降に、ラベル（文字型）、数値（200列）の組が、300000行　並んでいます。

出力したいものは、
１行目に記載されている表のサイズ（300000 200）をもつ、２行目以降のデータのmatrixです。
・ラベル（文字型）を、rownamesに格納
・colnames は、paste0("V",1:200)
・数値を数値型で、200列

２行目以降の表のサイズが大きいので、効率的な処理もできればと思います。

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

Re: Ｒ　レイアウトが異なるテキストファイルの読み込み
　　投稿者：青木繁伸 2017/05/26(Fri) 18:27 No. 22361
1 行目を読み飛ばせば，あとは，1行目に書いてあるサイズ情報を利用せずとも，read.table がちゃんと読んでくれます。
read.table(ファイル名, skip=1)
colnames もデフォルトで付けてくれます。
Re: Ｒ　レイアウトが異なるテキストファイルの読み込み
　　投稿者：明石 2017/05/26(Fri) 18:36 No. 22362
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

ご教示をいただきまして、大変に助かりました。
ありがとうございました。
Re: Ｒ　レイアウトが異なるテキストファイルの読み込み
　　投稿者：明石 2017/05/27(Sat) 15:23 No. 22366
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

昨日はご教示をいただき、大変に助かりました。
ありがとうございました。

ーーー

昨日の画像ファイルのデータ例について、追加のご質問がございます。　
どうぞよろしくお願いいたします。

1行目の　300000 200　は、２行目以降の表データのサイズ（行　列）を示しています。

「文字」＋「数値（200列）」の組が、300000行の繰り返しがあることを示しています。

（表のサイズは、題材により異なります。　固定ではありません）

２行目以降の表データのサイズ（行　列）を、１行目で明示している理由は、
ケースによっては、数値（200列）の後半に欠損がある場合を想定してからです。

件名の「レイアウトが異なるテキストファイルの読み込み」の真意は、ここにあります。

先生にご教示いただきたいことは、以下です。
お手数をおかけいたします。

・２行目以降に読み込む表データを、matrix型で格納する

・読み込んだ「文字」は、rownamesに格納する

・読み込んだ「数値（200列）」は、1列〜200列に数値データとして格納する。
　
・ケースによっては、数値（200列）の後半に欠損がある場合への対応として、
　データに欠損がある行については、列数が合うように（この場合は 200）
　その箇所をNA で補う。

お手数をおかけいたします。

ご教示いただければ、大変に助かります。

どうぞ、よろしくお願いいたします。
Re: Ｒ　レイアウトが異なるテキストファイルの読み込み
　　投稿者：青木繁伸 2017/05/27(Sat) 16:05 No. 22367
何度かお願いしておりましたが，添付図ではなく，仕様の全てを含む簡単なテストデータと，期待する結果をそえて質問してください。

説明文は曖昧さを含みますし，読み取りが不十分なこともありますから。
その点，データと結果があれば，仕様を満たすプログラムは一意に定まりますから。

「ケースによっては、数値（200列）の後半に欠損がある」というのは，一行ごとに数値の個数が違うことがあるということですか？
以下のようなものでよいですか？
file: "test.dat"
＝＝＝＝＝次行から＝＝＝＝＝
4 6
a 1 2 3 4 5
b 11 12 13
c 21 22 23 24 25 26
d 31 32 33 34
＝＝＝＝＝前行まで＝＝＝＝＝

func = function(fileName) {
	s = readLines(fileName, 1)
	s = as.numeric(unlist(strsplit(s, " ")))
	ncol = s[2]
	d = read.table(fileName, skip=1, row.names=1, header=FALSE, col.names=c("name", paste0("V", 1:ncol)), fill= TRUE)
	as.matrix(d)
}
func("test.dat")

結果
> a = func("test.dat")
> a
  V1 V2 V3 V4 V5 V6
a  1  2  3  4  5 NA
b 11 12 13 NA NA NA
c 21 22 23 24 25 26
d 31 32 33 34 NA NA
> class(a)
[1] "matrix"
【御礼】　Re: Ｒ　レイアウトが異なるテキストファイルの読み込み
　　投稿者：明石 2017/05/27(Sat) 17:00 No. 22368
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

今回もご丁寧にご教示いただき、誠にありがとうございます。
御礼を申し上げます。

先生からご指摘を受けましたように、
今後、質問をさせていただく際には、お書きになられたことに沿うようにいたします。

不愉快な思いをさせてしまい、大変に申し訳ございませんでした。
深くお詫びをいたします。

ーーー

今回の質問の内容については、先生がお書きになられた通りです。
私が期待していました通りです。

青木先生からのご教示は、大変に良いお手本となります。
さっそく、勉強をさせていただきます。

今回もありがたいご教示に感謝いたします。
ありがとうございました。


連続変数の自由度
　　投稿者：統計初心者 2017/05/17(Wed) 17:56 No. 22352
お世話になります。

一般線形モデルで、従属変数に連続変数のパラメーターを入れて検定したとき（サンプル数は48）、連続変数のパラメーターの自由度が１になるのはなぜでしょうか？
サンプル数は48なので、48-1=47になると思うのですが。
統計ソフトはJMP13を使用しています。

お忙しいところ恐れ入りますが、ご回答よろしくお願いいたします。
Re: 連続変数の自由度
　　投稿者：統計初心者 2017/05/18(Thu) 15:56 No. 22353
画像を追加します。
図の赤丸で囲った部分です。

Re: 連続変数の自由度
　　投稿者：青木繁伸 2017/05/19(Fri) 08:07 No. 22354
この類の質問は，回答に困りますね。
変数の不偏分散の自由度（n-1）と勘違いされているのでは？
Re: 連続変数の自由度
　　投稿者：統計初心者 2017/05/19(Fri) 11:46 No. 22355
青木先生

ご返答ありがとうございます。
この場合だと、線形モデルの式が y = xβ + e で、xが決まればyも決まるから、パラメーターの自由度は１という理解で大丈夫ですか？
Re: 連続変数の自由度
　　投稿者：青木繁伸 2017/05/19(Fri) 13:19 No. 22356
大元から掘り起こして理解するということもたいせつとは思いますが，「そういう定義である」と割り切って先に進む方が建設的では？

http://aoki2.si.gunma-u.ac.jp/lecture/Regression/mreg/mreg4.html
にある分散分析表の構成法と解釈法を理解することのほうがたいせつでは？
なお，そこの表ではp個の独立変数を持つ重回帰分析なので，一変量の単回帰分析なら p=1 です。重回帰分析の場合も，個々の独立変数の自由度は1です。
Re: 連続変数の自由度
　　投稿者：統計初心者 2017/05/19(Fri) 18:37 No. 22357
青木先生

大変参考になりました。ありがとうございます。


Ｒ　ポアソン回帰分析　glm関数の引数offsetの使い方
　　投稿者：明石 2017/05/08(Mon) 22:03 No. 22349
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

先日は、ポアソン分布のパラメータ計算方法について、分かりやすくご教示をいただき、誠にありがとうございました。

この投稿を契機に、Ｒのglm関数を用いて、ポアソン回帰分析の勉強を進めています。

glm関数のフォーミューラ（目的変数~説明変数）の記述方法は理解できました。

しかし、引数offsetの使い方が分かりません。

説明変数のうち、
どの変数をフォーミューラの右辺に書き、
どの変数を引数offsetに書くのか、
その使い分けが分かりません。

分かりやすい例示をしていただけると、助かります。
お手数をおかけいたします。
どうぞ、よろしくお願いいたします。
Re: Ｒ　ポアソン回帰分析　glm関数の引数offsetの使い方
　　投稿者：青木繁伸 2017/05/09(Tue) 06:20 No. 22350
example(glm) の中程以降に，family = gaussian ではありますが，offset を使う例が出てきますが...
glm> anorex.1 <- glm(Postwt ~ Prewt + Treat + offset(Prewt),
glm+                 family = gaussian, data = anorexia)
"R glm poisson offset" の 4 語で検索すると，いくつか出てきますが，
講義のーと: データ解析のための統計モデリング
eprints.lib.hokudai.ac.jp/dspace/bitstream/2115/49477/5/kubostat2008d.pdf
久保拓弥 著 - ‎2008
などにも例があるようです。
Re: Ｒ　ポアソン回帰分析　glm関数の引数offsetの使い方
　　投稿者：明石 2017/05/09(Tue) 09:05 No. 22351
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

ご紹介してくださいました資料を参照します。
お手数をお掛けして、申し訳ございませんでした。

いつも有難いご教示をいただき、大変に感謝しています。


共分散構造分析の重決定係数につきまして
　　投稿者：YUYU 2017/05/04(Thu) 20:04 No. 22341
青木先生

こんばんは。
はじめて質問させていただきます。YUYUと申します。

子育て中の看護師のバーンアウト因果モデルを作成しています。
ある集団に、この因果モデルを当てはめたところ、標準化推定値の重決定係数が、0.97となりました。適合度は、GFI = .925, AGFI = .890, RMSEA = .049でした。

ほかの集団での重決定係数は、おおむね、0.60から0.70でした。
0.97は、高すぎると思うのですが、いかがでしょうか。
ぜひ、先生のご意見をお伺いいたしたく、メールをさせていただきました。

どうぞよろしくお願いいたします。
Re: 共分散構造分析の重決定係数につきまして
　　投稿者：青木繁伸 2017/05/05(Fri) 08:39 No. 22342
0.97 という一つの数値だけ示されても，何ともいえないですね
Re: 共分散構造分析の重決定係数につきまして
　　投稿者：YUYU 2017/05/06(Sat) 10:04 No. 22343
青木先生

おはようございます。
お返事をいただきましてどうもありがとうございます。

当該の因果モデルを添付させていただこうと思ったのですが、どうやってもアップロードできませんでした。

一点、質問させてください。相関係数が高すぎると同質のものであることが疑われるといわれるように、重決定係数にも、一般的におかしいと考えられる値があるのでしょうか？
Re: 共分散構造分析の重決定係数につきまして
　　投稿者：青木繁伸 2017/05/06(Sat) 11:54 No. 22345
> どうやってもアップロードできませんでした。

本文中にテキストとして挿入するだけでよいのでは？

> 重決定係数にも、一般的におかしいと考えられる値があるのでしょうか？

何をもっておかしいとするか，に尽きるでしょう。

おかしな結果が出た理由は，分析結果ではなくデータそのものにある。
以下のような数値だけの結果を見せられても，おかしいかどうかわからないでしょう。
とても優れた予測式であるということかも知れない。
> x1 = c(3, 2, 1, 4, 2, 5, 6)
> x2 = c(4, 1, 2, 5, 4, 7, 6)
> y = c(17, 7, 9, 23, 15, 27, 35)
> a = lm(y ~ x1 + x2)
> summary(a)

Call:
lm(formula = y ~ x1 + x2)

Residuals:
       1        2        3        4        5        6        7 
-0.69528 -2.30901  1.96567 -0.02575  1.10730 -2.88412  2.84120 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)
(Intercept)   0.1760     2.3277   0.076   0.9434
x1            3.8026     1.1733   3.241   0.0316
x2            1.5279     0.9979   1.531   0.2005

Residual standard error: 2.612 on 4 degrees of freedom
Multiple R-squared:  0.9545,	Adjusted R-squared:  0.9318 
F-statistic: 41.96 on 2 and 4 DF,  p-value: 0.00207
以下に添付するのは，とてもよい（？）予測式が得られたというデータ例。標高の係数が絶妙。

Re: 共分散構造分析の重決定係数につきまして
　　投稿者：YUYU 2017/05/06(Sat) 17:28 No. 22346
青木先生

ご回答いただきまして、どうもありがとうございます。
大変勉強になります。

因果モデルのテキストへの貼り付けの件、試してみましたが張り付きませんでした。
図として添付してみます。

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

Re: 共分散構造分析の重決定係数につきまして
　　投稿者：青木繁伸 2017/05/06(Sat) 19:14 No. 22347
テキストの結果でなければテキストに挿入することはそもそも不可能ですが。

変数，潜在因子が多すぎるのではないかと思います。重回帰分析でも，独立変数をどんどん増やしていけば決定係数は大きくなりますので（自由度調整済み決定係数ならそのようなことはない）
解釈可能ならそれで良いかも知れませんが，より簡潔なモデルも検討してみたらいかがですか。
Re: 共分散構造分析の重決定係数につきまして
　　投稿者：YUYU 2017/05/07(Sun) 09:03 No. 22348
青木先生

おはようございます。
先生のアドバイスを参考に、他モデルも検討してみようと思います。

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


ポアソン分布の「プロイセン陸軍で馬に蹴られて死亡した兵士数」の事例
　　投稿者：明石 2017/05/03(Wed) 14:18 No. 22338
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

今回は、Ｒではなく、ポアソン分布について調べ物をしております。

Google先生に色々とお聞きしましたが分かりませんでしたので、投稿させていただきます。
どうぞよろしくお願いします。

−−−−−−−−−−
ポアソン分布の歴史的に有名な事例として、
「ロシアのプロイセン陸軍で馬に蹴られて死亡した兵士数」の事例があります。

Wikipedia他サイトでは、以下のように記載されています。
「プロイセン陸軍の14の騎兵連隊の中で、
　1875年から1894年にかけての20年間で馬に蹴られて死亡する兵士の数について
　調査しており、1年間当たりに換算した当該事案の発生件数の分布が
　パラメータ0.61のポアソン分布によく従うことを示している。」

これ以上に詳しい情報を得ることはできませんでした。

パラメータ0.61をどのように見つけたのか、計算したのか、適合度を確認したのか、
その方法に強い興味があります。

もし、これについて参考になる情報源をご存知でしたら、教えていただけましたらとても助かります。

お手数をおかけいたします。
どうぞ、よろしくお願いします。
失礼いたします。
Re: ポアソン分布の「プロイセン陸軍で馬に蹴られて死亡した兵士数」の事例
　　投稿者：青木繁伸 2017/05/03(Wed) 17:40 No. 22339
ポアソン分布のパラメータは「期待値（平均値）」です。
> 死亡した兵士の数 = 0:4
> 軍団数 = c(109, 65, 22, 3, 1)
> mean(rep(死亡した兵士の数, 軍団数))
[1] 0.61

> rbind(死亡した兵士の数, 軍団数)
                 [,1] [,2] [,3] [,4] [,5]
死亡した兵士の数    0    1    2    3    4
軍団数            109   65   22    3    1
> sum(軍団数)
[1] 200
> (0*109+1*65+2*22+3*3+4*1) / 200
[1] 0.61
Re: ポアソン分布の「プロイセン陸軍で馬に蹴られて死亡した兵士数」の事例
　　投稿者：明石 2017/05/03(Wed) 18:43 No. 22340
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

お休みのところ、ご丁寧にご教示をいただきまして、誠にありがとうございました。

一般化線形モデルのポアソン回帰モデルで、パラメータを推定するのかと思っておりました。

よい勉強になりました。
ありがとうございました。


尺度開発　テストー再テストの2回の集団の代表性を示す意味
　　投稿者：kashio 2017/04/28(Fri) 19:10 No. 22332
いつもご助言ありがとうございます

尺度開発の論文を投稿しました。査読者から以下４の指摘がありました。

4.研究結果
�‖仂櫃瞭胆�について、表1には1回目調査の結果のみですが、2回目調査に協力した集団（152名）は1回目調査の集団（774名）と特性が変わらないという前提でしょうか。1回目と2回目の調査対象に違いがないことを示す必要があります（サンプル代表性の担保）。つまり、表1には1回目と2回目の調査対象に関する記述集計を併記し、群間差があるか検定結果を示してください。検定方法については、p7のD.分析方法へ加筆してください。

これがどうもわかりません。

ほとんどの病院が２回（テスト―リテスト）はいやというので
1回目の調査の中の3病院が共同研究者の病院でここに頼みました。
300人ほどに頼みましたが、有効回答は半数以下でした。これもあり
テストーリテストは０・６で確かに低く課題が残りました。今後の課題だと書くつもりです。

2回目のテストは、テスト再テスト法にのみ使用します。
人数が少ないので、年齢、子供の有無、教育背景などにも差が出ますが、
示す意味が分からず、差があるから、なんだというのかもわかりません。
1回目の因子分析等の集団と再テストの集団と母集団が違うとだめ
違うのでしょうか？
ほかの尺度開発の論文をいくつか読みましたが、2回目の特性を示しているものは一つもありません。2回目は、1回目の一部です。そしてその同じ集団でテストーリテストの解析です。
代表性を示す意味も分かりません。
Re: 尺度開発　テストー再テストの2回の集団の代表性を示す意味
　　投稿者：青木繁伸 2017/04/29(Sat) 22:39 No. 22334
極端な場合を考えてみましょう。

尺度構成は８００人（男女４００人ずつ）のデータで行った。
２回目の調査は（極端にもほどがあるが）女１５０人に行った。
再現性の検討は２回とも調査できた女１５０人分のデータで行った。...
これでは，再現性があるなしをいえませんね。

そのあたりの担保が欲しいのでは？つまり，あなたの２回の調査で，対象者の背景因子に大きな違いがないと言って欲しいのかな。

しかし，再現性の検討は，尺度構成に使った集団とは違う集団に対して（妥当性もチェックできるし），もう少しサンプルサイズが大きい調査を行う必要があると思います。サンプルサイズが小さいと，「再現性がある」のか「再現性がないのだけど，それを検出できないだけ」なのかわかりませんからね。再現性を言うためにはなるべく少ないサンプルサイズで調査するのが...などと。
Re: 尺度開発　テストー再テストの2回の集団の代表性を示す意味
　　投稿者：kashio 2017/04/29(Sat) 23:36 No. 22335
ありがとうございました。男性が少なかったので女性だけを対象としてやり直しました。
青木先生の指摘通り、テスト再テストの人数が少ないうえ、共同研究者の大学の関連病院3つにしたため、特性も差が出ました。尺度の平均には差が出なかったのですが、年齢、実務職種、最終学歴に差が出ました。もうどうにも考察しようがないでしょうか？
探索的因子分析、確証的因子分析の結果は良好でした。クロンバックαは、0.69から0.84でした。テストリテストは、の相関係数は、0.6なので妥当性、信頼性はおおむね良好であった とは言えないでしょうか？
Re: 尺度開発　テストー再テストの2回の集団の代表性を示す意味
　　投稿者：青木繁伸 2017/04/30(Sun) 17:46 No. 22336
まさかの「男性が少なかったので女性だけを対象としてやり直しました。」ですか。

調査デザインとして問題が多いでしょうね。
Re: 尺度開発　テストー再テストの2回の集団の代表性を示す意味
　　投稿者：kashio 2017/05/01(Mon) 05:17 No. 22337
お返事ありがとうございます。対象が看護師なので最初から2000人以上対象にしないと男性看護師までは、無理かと危惧されたのですが、2回の調査の協力病院が得られず、尺度の平均にも性差が出たので、査読者に指摘されやり直しました。テストー再テストも0.6だったのでつらいですが、信頼性に課題が残ったが妥当性は確保されたとして査読者に返します。いつもありがとうございます。


対応のあるｔ検定？
　　投稿者：いちご 2017/04/23(Sun) 15:07 No. 22329
お世話になっております。
初歩的な質問だと思うのですが、教えていただければ幸いです。

４件法の尺度（６７項目）に対して、自分のことを回答した群（Ａ群）とその人のことをよく知っている人がその人のことを回答した群（Ｂ群）の対応データが３００ペア分あります。
このとき、A群の回答とＢ群の回答が等しいかそうでないかを検討するのに、どんな分析を行えばよいのでしょうか。
対応のあるｔ検定かなと思ったのですが、違うような気もします。

どなたか、ご教示いただけないでしょうか。
よろしくお願いいたします。
Re: 対応のあるｔ検定？
　　投稿者：課長 2017/04/24(Mon) 18:52 No. 22331
4件法ということであれば、各項目についてウィルコクソンの符号付き順位和検定のほうがふさわしいかと思いますが、６７項目あるということは６７回も検定を行うということになりそうですね。それはそれで問題かと思いますねぇ。多変量分散分析か、因子分析して各因子の合計点でウィルコクソンの符号付き順位和検定をするか、ですね。ｔ検定はあまりよろしくないかと思います。


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

二値データと連続変数の相関関係をみる方法をご教示お願い申し上げます。
Re: 二値データと連続変数の相関関係について
　　投稿者：課長 2017/04/22(Sat) 10:14 No. 22326
仰るとおり、１，０に意味は無いかと思います。計算上便利なだけでしょう。しかし、０，１を取るデータの数は意味があるのではないでしょうか？実際の計算では０や１を取る値の数、そして連続データの平均などの値を使うと思います。点双列相関係数やスピアマンの相関係数でいいんじゃないでしょうか？ちょっと自身がありませんので他の方の回答を期待します。
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つしか取りうる値がないものも間隔尺度という考えは自分には全くなかったので、
学習になりました。

「点双列相関係数」で検索してさらに勉強させて頂きます。
ありがとうございました。
Re: 二値データと連続変数の相関関係について
　　投稿者：課長 2017/04/23(Sun) 15:33 No. 22330
＞成瀬様
中途半端な回答で申し訳ありません。

＞青木先生
勉強になりました。今後とも宜しくお願い致します。


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因子すべて含めたモデルで論文投稿してみます。
この度は大変ありがとうございました。


患者さんに４年間介入した結果を分析する方法を教えてください
　　投稿者：分析初心者 2017/03/27(Mon) 18:07 No. 22309
初めて投稿いたします。よろしくお願いいたします。
患者さんに介入（入院前と入院中）し、質問紙調査によってデータを収集しました。介入の結果が年々患者さんの意識や態度に変化を及ぼしているのか、いないのかを知りたいと思っております。私たちの介入以外に様々な要因の影響があるとは重々承知しておりますが、検定をしたいと思っております。他に収集している変数は、年齢、性別、診療科、入院病棟などです。
目的変数は、「危険な行動の理解度（５件法）」や「入院中に履いていた履物の種類（３種類、質的データ）」など
意識や行動です。調査は、ある一時点の同じ時期・期間の患者さんを対象に調査しました。しかし、毎年患者さんは異なり、対応のあるデータではありません。介入前、介入１年後、２年後、３年後４年後とそれぞれ400名程度収集しております。このような場合、どのような分析方法が良いのがわかりません。分析方法についてお教えいただけますと幸いです。どうぞよろしくお願いいたします。
Re: 患者さんに４年間介入した結果を分析する方法を教えてください
　　投稿者：課長 2017/03/28(Tue) 12:16 No. 22310
「危険な行動の理解度」が介入前後で異なるかどうかなら、マンホイットニーのU検定になるんじゃないでしょうか？経年別に異なるかどうかならクラスカルワリス検定かなと思います。「履物の種類」が介入前後で異なるかどうかなら、カイ二乗検定ですかね。「危険な行動の理解度」や「履物の種類」を従属変数にして、介入や年齢、性別、診療科などのデータを独立変数とした回帰（重回帰や多重ロジスティックなど）も考えられます。いろいろ考えられますが、それぞれにいろいろ問題があるかと思います。
Re: 患者さんに４年間介入した結果を分析する方法を教えてください
　　投稿者：分析初心者 2017/03/28(Tue) 17:49 No. 22311
お忙しい中ご助言いただきましてありがとうございます。
＞経年別に異なるかどうかならクラスカルワリス検定かなと思います。
経年別にどう変化しているのか、その変化は有意なのかを調べたいと思っております。クラスカルワリス検定は、全く理解しておりませんので、勉強したいと思います。
＞いろいろ考えられますが、それぞれにいろいろ問題があるかと思います。
申し訳ありません。その通りだと思っております。
年齢が20歳から95歳までの幅がありましたので、中央値を基準に50才代以下と60才代以上に分けて分析しようと思います。
１年分のデータのみを分析対象としますと目的変数が量的なデータであれば重回帰分析など、「入院中に履いていた履物」のように質的なデータであれば、多重ロジスティック回帰分析ができるということだと理解いたしました。
青木先生の「統計手法の関連図」を順にしていきますと、「クラスカルワリス検定」に行き着きました。
検定の対象が「比率」も「平均値」どちらもありますので、あらゆるものを参考に分析したいと思います。ありがとうございました。
Re: 患者さんに４年間介入した結果を分析する方法を教えてください
　　投稿者：課長 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
情報ありがとうございます。


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

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

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

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

ーーー

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

タイタニックのデータをデータフレーム 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: Ｒ　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: Ｒ　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日

そこで以下についてお教え頂けないでしょうか。
ｔ検定を行う前提で、
サンプルサイズを計算するとき
予め、両群の日数の差と標準偏差が必要ということを理解しています。
この時の標準偏差を求める方法をお教え頂けないでしょうか。
恐らく、±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
青木先生
ご多忙の折、早々とご回答ありがとうございました。
とても勉強になりました。
とても感謝しております。
お礼まで


Ｒ　状態遷移を表形式で整理（２）
　　投稿者：明石 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"　のように、
遷移前のコンテンツ、遷移後のコンテンツを、"_"で連結した文字列を出力します。

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

Re: Ｒ　状態遷移を表形式で整理（２）
　　投稿者：青木繁伸 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: Ｒ　状態遷移を表形式で整理（２）
　　投稿者：明石 2017/03/06(Mon) 09:29 No. 22295
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

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

確認いたしました。

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

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


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

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

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

ーーー

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

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

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 n03,,,n07 の文字列は、コンテンツの名前です。

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

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

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

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

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

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

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

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

Re: Ｒ　状態遷移を表形式で整理
　　投稿者：青木繁伸 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: Ｒ　状態遷移を表形式で整理
　　投稿者：明石 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: Ｒ　状態遷移を表形式で整理
　　投稿者：青木繁伸 2017/03/05(Sun) 09:39 No. 22291
4 行目がおかしかったようです（関数名 names が落ちてしまった)
nameOfContents <- names(table(unlist(s)))
Re: Ｒ　状態遷移を表形式で整理
　　投稿者：明石 2017/03/05(Sun) 10:13 No. 22292
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

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

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

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

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

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


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

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

私は現在、ある上皮癌を浸潤群と非浸潤群の２群間について治療法や診断法についての検討を行っています。
この上皮癌には複数の組織型（１７種類）が含まれます。
組織型　浸潤群　非浸潤群
　
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
２群間における組織型の症例数の割合に差があるかどうかを調べたいので、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
青木先生、
いつもお世話になり、ありがとうございます、明石と申します。

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

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

ーーー

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

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

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

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

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

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

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

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

どうぞよろしくお願いいたします。
失礼いたします。
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が３種類の薬剤処理、REPが反復、Disは死滅株数、Healthは健全株数の１因子乱塊法実験です。
一般化線形モデルを用いて以下のように、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
　レスがつかないので、やはり難しいのですね。
　あまりスマートな方法ではないのですが、２つの処理ごとに以下のような分析を行うことにしました。
　普通の線形モデルで、２つの処理ごとにｔ検定を繰り返し適用するような感じですね。

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でＯＫのように読みましたが、
さらに、期待値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に５未満の桝目が全体の桝目の数の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
太郎様

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


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

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

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

ーーー

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

#----------------------------------------
#　年齢　⇒　年代（１０階級）
#----------------------------------------
#　　-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-"
#----------------------------------------

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

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 で指定した区間にデータがない場合には、エラーになるようです。

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

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

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

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

ご教示いただければ、大変に助かります。
何卒どうぞ、よろしくお願いいたします。
Re: Ｒ　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: Ｒ　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
青木先生

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

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

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

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