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

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


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

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


VBAでの購買間隔の平均を計算
  投稿者:波音 2019/09/20(Fri) 12:27 No. 22839
業務上の事情があってExcel VBAで「購買月間隔の平均」を計算するユーザー定義関数を作っています(作りました)。
自分でチェックしたところ、正しく動作はしているようですが。。。
記憶をたよりに書いていることもあり、何かもう少しスマートに書ける気がするようなしないような。。。

という状況ですので、もう少しスマートに書く方法があればご教示いただきたいです。

やりたいことは以下のようなものです。

例えば、私が毎月購買していれば4月〜3月までの12セルすべてに「1」が入ります。
もし購買していない月があれば「空白」が入ります。
仮に4,5,6;9,10;3に購買しなかった場合、空白は3, 2, 1ということになるので、その平均は(3 + 2 + 1) / 3という計算です。

この計算結果は2.00となりますが、意味としては「この人は平均的に2ヶ月間隔で買っていない(2ヶ月間隔でブランクがある)ということ。

VBAコード

Option Explicit

' 指定したセル範囲において、入力値と入力値の間の空白をカウントし、
'  空白数の平均を計算して返す関数

' ID×日付の購買頻度クロス集計表で、IDごとの購買間隔の日数の平均を求めるために作った関数

Function betweenAve(datRng As Range)

Dim dat As Variant ' 関数の引数で取得したアドレスにある数値を格納する配列
Dim n As Integer, i As Integer ' 配列の長さn、forの繰り返し用i

Dim blankCount As Integer ' 空白セルをカウントして格納していくハコ
Dim sumN As Integer ' 購買日間隔がいくつかるかカウントして格納していくハコ

Dim mySum As Double, ans As Double ' 空白セルを累計して格納するハコ

' まずは引数で指定されたセル範囲のアドレスから、セルに入力されている値をdatへ格納
' 二次元配列として格納されるので、行ベクトルの長さを取得
dat = Range(datRng.Address).Value
n = UBound(dat, 2)

' 各種変数を初期化しておく
blankCount = 0
sumN = 0
mySum = 0

' セルの値を1つずつチェックしていき、空白があればblankCountに「1」を格納していく
'  それをmySumへ累計していく
For i = 1 To n
If dat(1, i) = "" Then ' i番目のセルが空白がどうかチェック
blankCount = blankCount + 1 ' 空白なら「1」を代入
Else
If i <> 1 Then ' 空白でない場合、そのセルが最初のセルかどうかチェック
If dat(1, i - 1) = 0 Then ' 最初のセルでなければ、1つ前のセルをチェックして空白なら、
sumN = sumN + 1 '  sumNに+1カウント
End If
Else
'何も処理しない
End If

mySum = mySum + blankCount ' mySumへ空白の数を累計していく
blankCount = 0 ' 1ループの処理が終わるたびにblankCountを初期化する
End If
Next i

' 最終セルが空白の場合、ループ内で数がカウントされないのでループ後にカウントする
mySum = mySum + blankCount

' 最終セルが空白の場合、購買日の間隔数もカウントされないのでループ後にカウントする
If dat(1, n) = "" Then
sumN = sumN + 1
End If

' 購買日数間隔の平均を計算
ans = mySum / sumN

' 戻り値
betweenAve = ans

End Function


Re: VBAでの購買間隔の平均を計算
  投稿者:波音 2019/09/20(Fri) 12:29 No. 22840
Excelのワークシートのスクリーンショットも添付します。
このようなデータイメージです。


Re: VBAでの購買間隔の平均を計算
  投稿者:ts 2020/01/10(Fri) 18:02 No. 22883
Function 空白期間の平均(対象範囲 As Range)

Dim 前, 今, 分母

分母 = 0

前 = "最初なので、前はありません"

'--------------------------------

For Each 今 In 対象範囲

If 今 = "" And 前 <> "" Then

'空白期間の始まりだから

分母 = 分母 + 1

End If

前 = 今

Next
'--------------------------------

空白間隔の平均 = WorksheetFunction.CountBlank(対象範囲) / 分母

End Function

Re: VBAでの購買間隔の平均を計算
  投稿者:aoki 2020/01/10(Fri) 19:08 No. 22884
R だと,rle() を使えば簡単です。
data = c(
"111111111111",
"000110011110",
"110011001100",
"111111000000",
"101010101010")

for (str in data) {
str.vec = unlist(strsplit(str, ""))
result = rle(str.vec)
result = result$length[result$values == 0]
if (length(result) == 0) {
print(0)
} else {
print(mean(result))
}
}
結果
[1] 0
[1] 2
[1] 2
[1] 6
[1] 1
rle もソースを見れば簡単な関数です。余計なものを除けば以下のごとし。
rle2 = function(x) {
x = unlist(strsplit(str, ""))
y = x[-1] != x[-n]
i = c(which(y), n)
lengths = diff(c(0, i))
values = x[i]
return(list(lengths=lengths, values=values))
}
両方の言語を知っていれば,移植もそう難しくはないのでは?


R 区間分割処理
  投稿者:明石 2020/01/04(Sat) 19:15 No. 22879
青木先生 様;

新年あけまして、おめでとうございます。

お忙しいところを失礼いたします、明石と申します。
毎々、ご丁寧なご教示をいただき、誠にありがとうございます。

新年早々ですが、青木先生にご教示いただきたいことが出てきました。
年始のお忙しいところ、お手数をおかけいたします。
何卒どうぞよろしくお願いいたします。

----------------------------------------
AIC による,ヒストグラム(度数分布表)の最適階級分割の探索
http://aoki2.si.gunma-u.ac.jp/R/AIC-Histogram.html
の入力データを用いて、ご説明します。
x <- c(28.67, 40.29, 10.61, 33.85, 36.19, 20.63, 9.64, 15.26,
15.53, 73.62, 63.29, 32.77, 32.28, 11.90, 54.16, 4.73, 24.67,
17.66, 25.84, 22.89, 15.68, 5.48, 36.41, 20.33, 44.58, 57.23,
65.89, 57.91, 2.39, 9.15, 10.27, 3.04, 12.35, 32.78, 44.23,
31.14, 6.03, 27.90, 28.73, 42.09, 3.99, 9.74, 6.85, 0.16, 9.26,
7.72, 34.42, 32.77, 6.80, 10.45, 29.80, 5.89, 13.56, 50.55, 0.51,
0.19, 7.19, 5.94, 11.24, 32.32, 15.27, 29.64, 10.03, 2.01, 13.89,
20.83, 27.49, 14.46, 8.22, 27.81, 33.65, 38.57, 8.66, 1.40,
23.97, 15.11, 63.32, 7.76, 1.58, 48.66, 44.46, 0.02, 38.12,
18.51, 101.75, 34.16, 27.99, 5.22, 1.82, 8.22, 4.89, 97.50, 2.10,
26.19, 10.11, 8.39, 25.83, 1.05, 25.63, 18.35)

最適階級の区間幅が、以下のようにマトリクスで得られました。
mx
[,1] [,2]
[1,] -0.03000 10.68895
[2,] 10.68895 37.48632
[3,] 37.48632 64.28368
[4,] 64.28368 101.80000

findInterval(x, c(-0.03000, 10.68895, 37.48632, 64.28368))により、
入力ベクトルデータの各要素が所属する階級番号を得ることができました。
2 3 1 2 2 2 1 2 2 4 3 2 2 2 3 1 2 2 2 2 2 1 2 2 3 3 4 3 1 1 1 1 2 2
3 2 1 2 2 3 1 1 1 1 1 1 2 2 1 1 2 1 2 3 1 1 1 1 2 2 2 2 1 1 2 2 2 2
1 2 2 3 1 1 2 2 3 1 1 3 3 1 3 2 4 2 2 1 1 1 1 4 1 2 1 1 2 1 2 2

【ご教示いただきたい内容】
mxを人間が読み取り、
findInterval(x, c(-0.03000, 10.68895, 37.48632, 64.28368))
をコーディングするのではなく、
x、mxの2つを引数として渡すだけで、
findInterval(x, c(-0.03000, 10.68895, 37.48632, 64.28368))
を自動的に組み立てて、実行して、所属する階級を得ることができたら、大変に助かります。

最適階級化したいデータが多くあるので、
人手を介さないで自動処理ができたらと調べています。

eval()関数が利用できそうな気がして調べていますが、文字列の組み立てができません。

ご教示をいただけましたら大変に助かります。
新年早々、お手数をおかけいたします。
//

Re: R 区間分割処理
  投稿者:aoki 2020/01/05(Sun) 12:34 No. 22880
findInterval(x, mx[,1]) だけでよいのでは?

Re: R 区間分割処理
  投稿者:明石 2020/01/05(Sun) 16:26 No. 22881
青木先生 様;

新年あけまして、おめでとうございます。

お忙しいところを失礼いたします、明石と申します。

新年早々、ありがたいご教示をいただき、誠にありがとうございました。
今年もいい年になりそうな予感がしています。

eval()関数を見つけて、固執して、視野が狭くなっておりました。
//


「重複を除いた数」を意味する熟語
  投稿者:ts 2019/12/27(Fri) 15:07 No. 22869
統計学かどうか微妙な質問なのですが、

「重複を除いた数」を意味する熟語をご存知ないでしょうか。

「延べ人数」に対する「実人数」のような概念を示したいのですが、人数とは限りません。

「実数」と表現すると、自然数に対する実数(real number) と混同してしまいます。

「ユニーク数」と表現する人もいますが、ユニークという言葉は
「個性的な」という印象が先行するらしく、一般の人にうまく伝わりません。

うまい言い換えを知っていましたら、ぜひ教えてください。

Re: 「重複を除いた数」を意味する熟語
  投稿者:aoki 2019/12/28(Sat) 14:07 No. 22872
unique はいろいろな訳語がありますが,そのうちの一つとして「一意」があります。
「重複を除いた数」は「一意な要素数」,そのような要素は「一意な要素」でしょう。

Re: 「重複を除いた数」を意味する熟語
  投稿者:ts 2019/12/30(Mon) 15:50 No. 22878
ありがとうございます! 
「一意」なら、他の意味に解釈しようがないので、使いやすいです。
使わせていただきます。


重回帰分析の結果解釈について
  投稿者:山本 2019/12/28(Sat) 16:18 No. 22873
青木先生

はじめまして。
山本と申します。
統計の初心者で初歩的なご質問で恐縮でございますが、宜しくお願い申し上げます。

がん患者において、服用薬剤数が患者のADLに与える影響を調査しています。
当方の仮説は、入院中の薬剤数の減少は退院時のADLにプラスに働くと考えています。
入院中に薬剤数が減少した群を減少群、維持・増加した群を増加群に対象患者を分類し、単変量解析を行いました。
目的変数を退院時ADL(Bathel Index)、単変量解析の結果から有意因子を説明変数とした重回帰分析を行ったところ、服用薬剤数の変化量と退院時ADLに有意な関係が認められました[β(95%CI):-1.26(-2.21〜-0.32)]。
服用薬剤数の変化量は、退院時の薬剤数から入院時薬剤数を引いた値です。

この結果の解釈として、薬剤数が減少すると退院時ADLにプラスに働くと言えるのでしょうか?
御多忙のところ、申し訳ございませんが、ご教示の程、宜しくお願い申し上げます。

Re: 重回帰分析の結果解釈について
  投稿者:aoki 2019/12/28(Sat) 17:20 No. 22874
> 薬剤数が減少すると退院時ADLにプラスに働くと言えるのでしょうか?

あなたは,どう思うのですか?
そして,それのどこが不安なんですか?

Re: 重回帰分析の結果解釈について
  投稿者:山本 2019/12/29(Sun) 16:59 No. 22876
青木先生

ご返信ありがとうございました。
個人的には、薬剤数の変化量とADLとに負の関連性が認められておりますので、
自分の仮説通りの結果となっているのではないかと考えています。
しかし、自分の結果の解釈が正しいのかを相談できる人がいないため、ご相談させて頂きました。
宜しくお願い申し上げます。

Re: 重回帰分析の結果解釈について
  投稿者:aoki 2019/12/29(Sun) 22:33 No. 22877
> 自分の結果の解釈が正しいのか

分析法の実際のデータでの分析例,およびその分析結果の記述例,さらにはそれを論文という観点で述べた例も数多くあると思いますので,それらを参照されれば良いかと思います。
ネット上の話ではなく,古典的な手法ではあるが学会誌の参照など,同じ分野の同じテーマでの投稿例は参考になる,ということではなく「参考にすべき」ことだと思います。


R 時系列データの変換
  投稿者:明石 2019/12/28(Sat) 12:15 No. 22870
青木先生 様;

お忙しいところを失礼いたします、明石と申します。
毎々、ご丁寧なご教示をいただき、誠にありがとうございます。
改めて、御礼を申し上げます。

青木先生にご教示いただきたいことが出てきました。
年末のお忙しいところ、お手数をおかけいたします。
何卒どうぞよろしくお願いいたします。

----------------------------------------
時系列データで、Rに組み込みのAirPassengersがあります。
以下のようになります。

> class(AirPassengers)
[1] "ts"

> AirPassengers
Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec
1949 112 118 132 129 121 135 148 148 136 119 104 118
1950 115 126 141 135 125 149 170 170 158 133 114 140
1951 145 150 178 163 172 178 199 199 184 162 146 166
1952 171 180 193 181 183 218 230 242 209 191 172 194
1953 196 196 236 235 229 243 264 272 237 211 180 201
1954 204 188 235 227 234 264 302 293 259 229 203 229
1955 242 233 267 269 270 315 364 347 312 274 237 278
1956 284 277 317 313 318 374 413 405 355 306 271 306
1957 315 301 356 348 355 422 465 467 404 347 305 336
1958 340 318 362 348 363 435 491 505 404 359 310 337
1959 360 342 406 396 420 472 548 559 463 407 362 405
1960 417 391 419 461 472 535 622 606 508 461 390 432

これを、単純にデータフレーム形式に変換すると、以下のようになります。
日付データが消失してしまいます。

> head( as.data.frame(AirPassengers) )
x
1 112
2 118
3 132
4 129
5 121
6 135

例えば、以下のイメージ例のように、
日付データを保持したまま、
データフレームに変換できる方法があれば、
ご教示をいただけましたら、大変に助かります。

(出来上がりイメージ例)

date AirPassengers
1 1948-02-01 112
2 1949-02-01 118
3 1949-03-01 132
4 1949-04-01 129
5 1949-05-01 121
6 1949-06-01 135


年末のお忙しいところ、お手数をおかけいたします。
//明石

Re: R 時系列データの変換
  投稿者:aoki 2019/12/28(Sat) 13:39 No. 22871
何かあるはずと,seq のヘルプを眺めて,以下のようにすれば良いとわかりました。
> data.frame(date = seq(as.Date("1949/01/01"), as.Date("1960/12/01"), by="month"), AirPassengers = AirPassengers)
date AirPassengers
1 1949-01-01 112
2 1949-02-01 118
3 1949-03-01 132
4 1949-04-01 129
:
139 1960-07-01 622
140 1960-08-01 606
141 1960-09-01 508
142 1960-10-01 461
143 1960-11-01 390
144 1960-12-01 432

御礼(Re: R 時系列データの変換)
  投稿者:明石 2019/12/28(Sat) 18:26 No. 22875
青木先生 様;

お忙しいところを失礼いたします、明石と申します。

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

Rコードだけでなく、着眼点は大変に良い勉強となりました。

拡張パッケージに頼らなくても、本体に、ちゃんとあるのですね。

おかげさまで、良い年を迎えることができます。
//


サーストンの一対比較法で求めた間隔尺度から有意差の有無を調べる方法
  投稿者:Nakamori 2019/12/26(Thu) 02:46 No. 22866
青木先生

初めまして,大学生の中森と申します.
表題の件ですが,どのように検定をすればよいか分からないため,
ご教示いただけたらと思います.

--以下事例--
現在制御手法の有効性について調べており,硬さの異なる5種類のばねを
4つの制御手法も用いて,どの制御手法が最も硬さを判別できるか検証したいと考えております.

サーストンの一対比較法を用いて,制御手法1つずつばねの硬さの間隔尺度を算出し,
どの制御手法が一番理想的(100%の判別率)に近いかまでは算出出来ました(ユーザn=7).

しかし,制御手法に対してどのようにすれば有意差が示せるかが分かりません.

お忙しいところ恐縮ですが,ご教授頂けたら幸いです.

Re: サーストンの一対比較法で求めた間隔尺度から有意差の有無を調べる方法
  投稿者:aoki 2019/12/26(Thu) 13:59 No. 22867
サーストンの一対比較法は検定手法ではないので,検定は出来ません。
実際にどのようなデータを取ったのか今ひとつはっきりしないのですが,二元配置分散分析とか乱塊法が適用できるか検討してみればいいかもしれません。

Re: サーストンの一対比較法で求めた間隔尺度から有意差の有無を調べる方法
  投稿者:Nakamori 2019/12/26(Thu) 20:20 No. 22868
青木先生

ご返信ありがとうございます.
要領を得ない質問で申し訳ありませんでした.
先生のおっしゃる通り,二元配置分散分析や乱塊法を検討してみます.

その時はまたお力になっていただけると幸いです.
重ねてお礼申し上げます.


AIC による,ヒストグラム(度数分布表)の最適階級分割の探索
  投稿者:明石 2019/12/09(Mon) 07:36 No. 22863
青木先生 様;

お忙しいところを失礼いたします、明石と申します。
毎々、ご丁寧なご教示をいただき、誠にありがとうございます。
改めて、御礼を申し上げます。

青木先生にご教示いただきたいことが出てきました。
何卒どうぞよろしくお願いいたします。

----------------------------------------

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

について、お聞きしたいことがございます。

> floor(2*sqrt(length(x))-1)
[1] 19

最初は、19階級です。

19階級が併合されて、
最適階級数は 4 となり、
その時のAICは、487.53 となります。

19 → 4 と階級数が大きく下がったので、
途中の階級数、ためしに、明示的に、階級数を7にして実行しました。

ans <- AIC.Histogram(x, 0.01, 7)

その時のAICは、297.84 となります。

こちらの方がAICが低いのが、気になります。

ans <- AIC.Histogram(x, 0.01, 7)などのような使い方は、間違った使い方でしょうか。

ご教示をいただければ助かります。
お手数をおかけいたします。

Re: AIC による,ヒストグラム(度数分布表)の最適階級分割の探索
  投稿者:aoki 2019/12/10(Tue) 01:19 No. 22864
最初の階級分けがあってのその後の探索ですので,出発点が異なれば異なった結果になる事は避けられないでしょう。
たとえば,極端な場合(サンプルサイズにもよりますが)出発点が100階級なんかに設定されていたら,そもそもそんな階級分けが不適切なら,両端の階級をまとめるという手順を踏むだけですから,結果も不適切になるのではないでしょうか?

Re: AIC による,ヒストグラム(度数分布表)の最適階級分割の探索
  投稿者:明石 2019/12/10(Tue) 09:36 No. 22865
青木先生 様;

お忙しいところを失礼いたします、明石と申します。

ご教示をいただき、ありがとうございました。
勉強になりました。

つまらない質問で申し訳ございませんでした。
今後、十分に気をつけます。

失礼いたします。


混合計画の分散分析についてです
  投稿者:初学者 2019/11/21(Thu) 13:20 No. 22861
大学の授業で 分散分析の解釈の課題が出ました。
薬物Aの投与による課題成績の影響を調べる実験です。
実験デザインは 対応なし(薬物A投与・投与なし)と対応あり(投与前・投与後)
の2要因混合計画です。
    投与前  投与後
投与群  3.3   3.6
統制群  2.4   2.4  という結果になり
群の主効果のみ5%水準で有意だとわかりました。
この場合の解釈についてなのですが
私としては 交互作用が有意で
投与群における時期の単純主効果と 投与後の群の単純主効果が有意なら
薬物Aには課題成績向上の効果の可能性あり と言えると思ったのですが

群の主効果のみだったので どう解釈すべきか分からなくなってしまいました
この主効果に 意味はあるのでしょうか?
また 投与前の得点に差があるのですが これは投与前の群の単純主効果が有意でなければ この差は偶然だった として 群分けに問題はなかったと考えてもよいのでしょうか?

先生のご教授のほどよろしくお願い致します。


データの対応の有無について
  投稿者:坂上 2019/11/05(Tue) 16:16 No. 22855
初めて投稿させていただきます。お手数をおかけしますが、ご教授いただけましたら幸いです。
質問紙調査において、

設問1
パンを毎日食べますか? はい、いいえの2択
設問2
設問1で「はい」の方のみお答えください。
食べているパンは自分で購入していますか? はい、いいえの2択
設問3
白米を毎日食べますか? はい、いいえの2択
設問4
設問3で「はい」の方のみお答えください。
食べている白米は自分で購入していますか? はい、いいえの2択

という設問(内容は一部改変しています)で、設問2と4で「はい」「いいえ」の割合が異なるかを調べたいと思っています。
この場合、設問2と4の回答者の中には、「設問2だけ答えた人」「設問4だけ答えた人」「設問2と4ともに答えた人」が混在していることになるかと思います。
この場合、設問2と4のデータは対応があるデータなのか、対応がないデータなのか(設問2と4ともに答えた人は対応があると思いますが)、どのように考えたらよいのでしょうか。

Re: データの対応の有無について
  投稿者:aoki 2019/11/05(Tue) 22:21 No. 22857
質問の条件が不備なため,あなたの知りたいことを得られる情報を十分に提示されてていません。

質問としては,「あなたはパンを毎日食べますか」,「あなたは白米を毎日食べますか」,「それは自費で購入していますか」の3項目を調べなければならないでしょう。もっとも,「パン」か,「白米」か「それ以外」かの3つの選択肢にしないといけないでしょうけど。「毎日」という強制用語も妥当かどうか疑問とか。

データを取った後で,ああだこうだいうのは,「無理です」。

- - - - - - - - -

話をわかりやすくするつもりで,例え話にして質問するということは質問する人にはよくあることですが,問われる方にしてみれば例え話に引きづられて本質から外れることはよくあることです。例え話が上手くないということ。あなたが思うほど,例え話でよく分かると言うことはありません。かえって,変な方向へ行ってしまうことが多いのです。

実質データはそのままで(分からない人には分かってもらわなくてもよい。だって,答えを出せる可能性は低いのだから。ごめん!),正確な背景をちゃんと述べることが正解に近づく近道ですよ(この原則が必須条件であることは,はもう何回となく繰り返し検証されている経験則なので)

Re: データの対応の有無について
  投稿者:坂上 2019/11/06(Wed) 10:58 No. 22858
青木先生

 稚拙な質問にも関わらずご返信をいただき、ありがとうございました。複数人で共同して行っている調査なので、個人の判断で内容をそのまま書いてよいものか、と思っている段階で投稿してしまいました。ご指摘のとおり、調査内容そのものでないとやりとりが的確でなくなってしまうこと、理解できました。とりあえず、今回の質問自体は取り下げさせていただきたく思います。不備のある質問にも関わらずご教授くださり、重ねて感謝申し上げます。


複数のヒストグラムの取り扱いについて
  投稿者:平沢 2019/09/20(Fri) 16:35 No. 22841
青木先生

初めての投稿になります。
大学の研究で行き詰ってしまったためご教授お願い致します。

試験体に含まれた薬剤量と、その薬剤の主成分が試験体のどこにどれほど入っているかを分析しています。データとしては、1つの試験体をおおよそ25000分割し(試験体毎に分割数は異なる)、その1点1点で主成分濃度を得ています。
例えば、試験体Aでは薬剤量149kg/m3、主成分濃度のサンプル数26510個といったものです。

以上のようなデータが30試験体分あり、それぞれで主成分濃度のヒストグラムを作成致しました。ヒストグラムの階級数はスタージェスの公式を用いて算出しています。
ヒストグラムの形について、主観では、薬剤量が 1)一定値以下 2)一定値付近 3)一定値以上の3タイプに分類されるように見えておりますが、これを客観的に言える方法が無く困っております。
2)一定値付近では比較的正規分布しているように見えたためコルモゴロフ=スミルノフ検定を行いましたが、帰無仮説は棄却されました。

稚拙な説明で申し訳ございませんが、ご助言いただければ幸いです。

Re: 複数のヒストグラムの取り扱いについて
  投稿者:青木繁伸 2019/09/20(Fri) 23:06 No. 22843
データのイメージがはっきりつかめないので,シミュレーションデータを作って分析してみました。

n 試験体(貴方の例では n = 30) のそれぞれについてかなりのデータがあり,それを class (貴方の例では class = 25000) 区間に分割してヒストグラム(度数分布)を作った。このヒストグラムについて,似ているもののグループ分けをしたい?
set.seed(123)
# class = 2500
class = 250 # ヒストグラムの階級数
N = 10000 # データ数
# n = 30
n = 15 # 試験体
data0 = matrix(0, n, N) # 大元のデータ(試験体数 x データ数)
data = matrix(0, n, class) # 試験体ごとの度数分布
Mean = runif(n, -1, 1) # 試験体ごとの平均値(取りあえず)
for (i in 1:n) { # 試験体ごとの元データ(取りあえず正規分布に従うとして)
data0[i, ] = rnorm(N, mean=Mean[i])
}
Min = min(data0) # データ全体の最小値(度数分布を作るために)
Max = max(data0) # データ全体の最大値(度数分布を作るために)
for (i in 1:n) { # 試験体ごとの度数分布(ヒストグラムは描かない)
data[i, ] = hist(data0[i,], breaks=seq(Min, Max, length=class+1), plot=FALSE)$counts
}
さて,このような,試験体 x 度数分布の集計結果をデータとみなして主成分分析を行うとどうなるか?
a = prcomp(data) # scale.=TRUE が必要かどうか?
主成分分析の結果を図示する。試験体の配置をみるために,主成分得点を図示する。
図をクリックすると拡大表示

次の発言に続く


Re: 複数のヒストグラムの取り扱いについて
  投稿者:青木繁伸 2019/09/20(Fri) 23:09 No. 22844
続き

第1主成分の順番で度数分布を図示してみる。
plot(c(Min, Max), c(1, n), type="n", xlab="class", ylab="sample no.", yaxt="n")
t = order(a$x[,1])
for (i in 1:30) {
points(seq(Min, Max, length=class+1), rep(i, class+1), pch = 19, cex=sqrt(data[t[i],] )/10)
text(-4, i, t[i], col="red")
}
先の図からも似ているもの同士の位置関係がわかるが,実際どのように似ているかを試験体ごとの度数分布を示してみる。
第一主成分の順番で表示してみると,平均値の大きさにしたがって配置されていることがわかる。


Re: 複数のヒストグラムの取り扱いについて
  投稿者:青木繁伸 2019/09/20(Fri) 23:15 No. 22846
第2主成分の順でプロットしてみる。
plot(c(Min, Max), c(1, n), type="n", xlab="class", ylab="sample no.", yaxt="n")
t = order(a$x[,2]) # 第2主成分を指定する
for (i in 1:30) {
points(seq(Min, Max, length=class+1), rep(i, class+1), pch = 19, cex=sqrt(data[t[i],] )/10)
text(-4, i, t[i], col="red")
}
右上,左上,下部の3パターンに分類できる。


Re: 複数のヒストグラムの取り扱いについて
  投稿者:青木繁伸 2019/09/20(Fri) 23:21 No. 22847
結局の所,試験体ごとの度数分布をデータと見なして似たもの同士を集める統計手法を使えばどうだろうかということ。

因子分析,対応分析,クラスター分析でもよいだろうと思う。(私が理解したデータ構造があたっていれば)

Re: 複数のヒストグラムの取り扱いについて
  投稿者:平沢 2019/09/22(Sun) 14:06 No. 22848
青木先生

ご返信ありがとうございます。
まだRを初めて間もないため理解に時間がかかり
返信が遅くなりました、申し訳ございません。

データ構造につきましては大丈夫です。
スケールはnが試験体名、classは薬剤濃度 kg/m^3のため
TRUEで良いかと思います。

まずは主成分分析をしてみる、というアドバイスありがとうございます。
素人には中々頭がそこまで回らず、先生にお聞きして本当に良かったです…。
実際のデータで試してみたいと思います。

またそこで何かありましたら頼りにさせていただければ幸いです。
ご多用のところ、ありがとうございました。

Re: 複数のヒストグラムの取り扱いについて
  投稿者:平沢 2019/09/24(Tue) 15:54 No. 22849
青木先生

度々失礼致します。
不明な点が2点ありましたのでお聞きしたく存じます。

>data0 = matrix(0, n, N) # 大元のデータ(試験体数 x データ数)
N数が試験体ごとに異なるのですが、
この場合どういう風に行列の指定をしたらよろしいでしょうか。
最もデータ数の多いものに合わせれば問題ないでしょうか。

>Mean = runif(n, -1, 1) # 試験体ごとの平均値(取りあえず)
なぜこの値が試験体毎の平均値であるのか、わかりませんでした。

以上、2点お願い致します。

Re: 複数のヒストグラムの取り扱いについて
  投稿者:青木繁伸 2019/09/25(Wed) 16:35 No. 22850
挙げたのは例に過ぎないので,実際の貴方のデータを使えば良いのです。

> N数が試験体ごとに異なるのですが、
> この場合どういう風に行列の指定をしたらよろしいでしょうか。

最終的な度数分布データを作るための架空データなので,貴方の場合は試験体ごとの実際のデータから度数分布データを求めればよいだけです。

> なぜこの値が試験体毎の平均値であるのか、わかりませんでした。

これも,最終的な度数分布データを作るための架空データなので,貴方の場合は試験体ごとの平均値はそれこそ試験体の特徴を表す実際の値があるでしょう。

ようするに,貴方のやることは,試験体ごとの度数分布データをまとめるだけです。以下のような 30 × 25000 の行列データ

階級1  階級2  …  階級25000
試験体1  ○○  ○○  …    ○○
試験体2  ○○  ○○  …    ○○
:
試験体30  ○○  ○○  …    ○○

Re: 複数のヒストグラムの取り扱いについて
  投稿者:平沢 2019/09/25(Wed) 18:33 No. 22852
青木先生

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

勘違いをしていました、納得できました。
失礼いたしました。

追加で質問したいのですが
複数のヒストグラムを1つにして比べたいとき、
適切な階級幅の設定方法で悩んでおります。

例えば、私が扱っているデータは
試験体A 最小値0、最大値7500、中央値2000、サンプル数27000
試験体B 最小値130、最大値10600、中央値2500、サンプル数25000
試験体C 最小値0、最大値7000、中央値1600、サンプル数26000   
etc…
といったオーダーです。

1つずつ作る際はスタージェスの公式を参考にした階級数から
最大値と最小値の差を割ってでた数値をキリの良い階級幅に設定しています。
しかし、例えば試験体BとCでは以上の方法で算出した階級幅には
大きく差が出てしまいます。
できる限り偏りのないグラフにするために
何か方法があればお願い致します。

※私用のため3日間は返信ができないこと、ご容赦ください。

Re: 複数のヒストグラムの取り扱いについて
  投稿者:青木繁伸 2019/09/25(Wed) 21:05 No. 22853
注意すべきは,全ての試験体で同じ階級設定をしないといけないこと。(当然ですが)

何回も書いているが,以下のようなデータ行列を用意すること
階級1  階級2  …  階級25000
試験体1  ○○  ○○  …    ○○
試験体2  ○○  ○○  …    ○○
:
試験体30  ○○  ○○  …    ○○

階級数は 25000 も必要ないとは思うが。

Re: 複数のヒストグラムの取り扱いについて
  投稿者:平沢 2019/09/30(Mon) 12:55 No. 22854
青木先生

階級幅の平均をとるなどして調整してみます。
何度もご回答いただきありがとうございました。


分散分析で大丈夫でしょうか
  投稿者:コロン 2019/09/18(Wed) 08:06 No. 22833
お世話になっております。ご教授ください。

3クラスに対して、授業前後の意識の変化を5件法で問うています。項目は5項目あり、それらの平均値を特性値にします。

この分析手法ですが、2元配置分散分析でよろしいでしょうか?それとも、五件法ということて、ノンパラの手法、フリードマンで行うべきでしょうか?

よろしくお願いいたします。

Re: 分散分析で大丈夫でしょうか
  投稿者:コロン 2019/09/18(Wed) 09:06 No. 22834
追加で書かせてください。

5項目の平均値を特性値にすると書きましたが、足し算した値を特性値にすることも考えましたが、大きな違いはないのでしょうか。

基本的な質問で申し訳ございません。

Re: 分散分析で大丈夫でしょうか
  投稿者:青木繁伸 2019/09/18(Wed) 19:05 No. 22835
> 平均値を特性値にすると書きましたが、足し算した値を特性値にすることも考えました

5 で割るか,そのままかで,scale の違いだけなのでどちらも同じ結果になるのでは?

Re: 分散分析で大丈夫でしょうか
  投稿者:コロン 2019/09/18(Wed) 19:18 No. 22836
青木先生

ありがとうございます。このデータは、どの手法を用いて差があるかどうかを見ればよろしいのでしょうか。

Re: 分散分析で大丈夫でしょうか
  投稿者:青木繁伸 2019/09/19(Thu) 08:14 No. 22837
5 件法でしたね。だとすると合計を取るのも平均値をとるのも不適切。
ということで,フリードマン検定か。

Re: 分散分析で大丈夫でしょうか
  投稿者:コロン 2019/09/19(Thu) 10:13 No. 22838
青木先生

お返事ありがとうございます。ということは、一つの概念を表していると考えている5項目を、それぞれにフリードマンという理解でよろしいでしょうか。


数量化粁爐砲弔まして
  投稿者:こだま 2019/09/11(Wed) 17:30 No. 22830
青木先生

お世話になります。
こだまと申します。

本掲示板の 2019/08/14(Wed) 13:13 No. 22795 にて
統計分析の手法(判別分析・数量化粁燹砲砲弔い銅遡笋気擦討い燭世い燭發里任后

先生からのご助言やその他の参考書等を読み、最終的には数量化粁爐砲栃析を進めた
のですが、ご教示いただきたい事があり、再度書き込みさせていただきます。

数量化粁爐砲董公園のあり/なしを目的変数として、その要因を分析したのですが、分析に使用したサンプル数が、あり:144、なし:5 でした。
相関比や判別的中率については、それぞれ0.49、96.0%となり、精度としても良さそうなので、この分析に基づいて文章を書こうかと思っていたのですが、今回のようにサンプル数に偏りがある場合に問題はないのだろうか(分析するデータとしてはよろしくない?)
と思いました。

ある書籍では、「数学的にはサンプルが3以上あればよい」というような記載していたのですが、本当にそれで良いのか、分析上問題がある場合はどのような分析ができそうかについてご教示いただけますと幸いです。

Re: 数量化粁爐砲弔まして
  投稿者:青木繁伸 2019/09/11(Wed) 21:08 No. 22831
「数学的にはサンプルが3以上あればよい」 という,「サンプル」が何を意味するか?そもそも,統計学で「サンプル」という用語を使うだけで,怪しさ満載ですが。

しっかりした分析プログラムで分析できれば,「数学的に問題なし」と理解してよいでしょう。「サンプル」がなにかを置いておいても,「サンプルの3個が全く同じ値を取ったら」数学的には問題ありというか,分析プログラムはエラーで停止するでしょう。

それはさておき,「数学的に問題なし」と「実質的に妥当性あり」とは違います。例えばあり144,なし0なら,数学的に問題あり。あり144,なし1なら「数学的に問題なし」ですが,「妥当性がある」とは言えないでしょう。あり144,なし2なら?あり144,なし3なら?...
どのあたりから妥当性があるといえるでしょうね?これは統計学で答えることはできない問題でしょう。大方の人が妥当だと思ってくれないと,論文としてアクセプトするのは難しいと言うことになるでしょう。

Re: 数量化粁爐砲弔まして
  投稿者:こだま 2019/09/11(Wed) 21:33 No. 22832
青木先生

お世話になっております。

御返事いただきありがとうございます。

ご指摘いただいたように、色々と当方の専門とする分野の論文をレビューしている中で、統計の面で本当に妥当なのか?というものがいくつかあり、妥当性というところでもやもやしておりました。

ここで質問してお答えいただく事ではなかったにも拘らず、ご教授いただきありがとうございます。大変勉強になりました。


教育心理 重回帰 ロジスティック
  投稿者:佐藤 2019/07/01(Mon) 16:11 No. 22770
教育心理分野で、ある実習後の学習意欲を目的変数とし、実習中の環境要因2つ、実習生の個人要因2つを説明変数として重回帰分析をしようと考えていました。

学習意欲に対して正規性の検定を行ったところ正規性が認められなかったため、重回帰分析を適応できませんでした。

1. この場合、ノンパラメトリックの相関分析までしか適応できないのでしょうか?

2. 学習意欲のデータを見たところ二峰性があったため、学習意欲の中央値で2群に分割し、高学習意欲群1と低学習意欲群0にダミー変数化してロジスティック回帰分析をしても良いのでしょうか?(それともこれはデータの恣意的な変更や情報量の削減という観点で許されないことなのでしょうか?)

統計に明るくなくおかしなことを聞いているかもしれませんが、ご教授頂けますと助かります。

Re: 教育心理 重回帰 ロジスティック
  投稿者:青木繁伸 2019/07/01(Mon) 21:10 No. 22771
> 学習意欲に対して正規性の検定を行ったところ正規性が認められなかったため、重回帰分析を適応できませんでした。

重回帰分析において,目的変数が正規分布するかどうかは関係ありません。というか,目的変数がどのような分布をするかによって,適切な分析方法があります。

例えば,一つ前の質問のような,目的変数が 0/1 の二値変数の場合にはロジスティック回帰分析とか,その他にも glm 関数は目的変数が様々な場合に対して分析を行うことができるように対処しています。
binomial(link = "logit")
gaussian(link = "identity")
Gamma(link = "inverse")
inverse.gaussian(link = "1/mu^2")
poisson(link = "log")
quasi(link = "identity", variance = "constant")
quasibinomial(link = "logit")
quasipoisson(link = "log")
また,単に重回帰分析の場合であっても,従属変数が正規分布に従うことを要件にしてはいません。

Re: 教育心理 重回帰 ロジスティック
  投稿者:佐藤 2019/07/02(Tue) 14:18 No. 22772
青木先生

ご返信有難うございます。
不勉強でした重回帰分析は従属変数の正規分布を前提としていないんですね。
ネット上でK-S検定を行って正規性を確認してから、
重回帰分析を行っている例があったので、重回帰分析は正規分布を仮定している統計法と誤解をしていました。

従属変数は間隔尺度のデータなので、無理に2値化してロジスティック回帰分析をするよりも、
重回帰分析を行なうほうがデータに沿った分析方法のように考えられます。

重回帰分析もう一度勉強し直してみます。

青木先生、今回は誠にありがとうございました。

Re: 教育心理 重回帰 ロジスティック
  投稿者:佐藤 2019/08/23(Fri) 11:48 No. 22812
青木先生

以前、丁寧にご回答頂きまことにありがとうございます。
以前の質問の続きなのですが、目的変数、説明変数ともに正規性の無いデータで重回帰分析を行なったところ、学会にて正規性の無いデータで重回帰分析は使用できないとの指摘を受けました。

「すくできる!リハビリテーション統計」という書籍や他医療統計の書籍に
「重回帰分析では独立変数も従属変数も正規分布した数量データ(間隔尺度か比率尺度)である必要があります。」

という記述があります。
以前のやり取りの中で、正規分布しないデータでも重回帰分析が適用できるという、
当方の理解が間違っていたのでしょうか?

Re: 教育心理 重回帰 ロジスティック
  投稿者:青木繁伸 2019/08/23(Fri) 14:21 No. 22813
その本の著者も,学会(査読者?)も間違えています。間違ったことを書いているクソページもたくさんありますね。

前に説明したことに基づいて,反論しなかったのですか?

説明変数は正規分布しなくても差し支えありません。ダミー変数を説明変数に使う場合を考えて見ればわかるでしょう。ダミー変数は 0/1 の二値データで,正規分布などしません。カテゴリーデータは重回帰分析に使えない(使うなら数量化 I 類を使いなさい)などという,とんでもない時代錯誤の考えを持っている人もまだいるのだろうか?また,実験計画によるデータの場合,説明変数は特定の値を取ります(薬剤投与量と効果の場合,薬剤投与量は例えば 10mg, 20mg, ..., 100mg などとされるでしょう。だれも薬剤投与量を正規乱数から選ぼうなどとは思わないでしょう)。

目的変数は,正規分布しなくても構いません。正規分布しないといけないのは残差の分布です(さらに,それさえも必須条件ではない)。
誤差が正規分布しないと,予測値の信頼区間を求めることができないからです。誤差(測定誤差も含む)は普通は正規分布に従うものです。なので,統計分析できるのです。

日本語のページは信頼できないということなら,英語のページを検索してはいかがでしょうか。

https://www.statisticssolutions.com/assumptions-of-multiple-linear-regression/
Multiple regression assumes that the residuals are normally distributed.

http://www.restore.ac.uk/srme/www/fac/soc/wie/research-new/srme/modules/mod3/3/index.html
Normally distributed residuals: The residuals should be normally distributed.

https://www.m3.com/open/clinical/news/article/604122/
回帰式の形は、被説明変数のタイプによって変えることができます。被説明変数が連続変量の場合は重回帰分析、2値変数(時間依存性がないイベントの発生または非発生、死亡または生存、治癒または非治癒など)の場合はロジスティック回帰分析、打ち切りのある2値変数の場合はCox回帰分析を用います。

重回帰分析では、残差(回帰分析による予測値と実測値の差分)の分布が正規分布である必要があります。

https://bellcurve.jp/statistics/course/9700.html
回帰モデルを考えるにあたって、誤差 μi にはいくつかの仮定している条件があります。
1. μi の期待値は0である
2. μiの分散は常にσ2 である
3. 異なる誤差 μi、μj は互いに独立である
これらの3条件から、μi は互いに独立に正規分布 N(0, σ2) に従うと仮定されます。

==================================

具体例でも示してみましょう。正規分布すべきは「残差」です。

理論モデル y = 2.5 + 4.8 * x

重回帰モデルは y = a + b * x + ε
εが正規分布

データを用意する

n = 1000
x = runif(n, min=0, max=100) # 独立変数
EPSILON = rnorm(n, mean=0, sd=10) # ε:正規分布
y = 2.5 + 4.8 *x + EPSILON # 誤差を含む従属変数

y は正規分布しない!!

hist(y, main="図1")

予測は十分

plot(y ~ x, pch=19, col="#aa000055", alpha=0.3, main="図2")
ans = lm(y ~ x)
summary(ans)
abline(ans)
text(10, 400, sprintf("y = a + b x\na = %.3f\nb = %.3f", ans$coef[1], ans$coef[2]), pos=4)

誤差(予測誤差)は正規分布する

plot(ans$residuals ~ x, pch=19, col="#aa000055", main="図3")

hist(ans$residuals, main="図4")

図は,クリックすることで拡大表示


Re: 教育心理 重回帰 ロジスティック
  投稿者:佐藤 2019/08/28(Wed) 09:41 No. 22827
青木先生

出先で返信書き込み遅くなり大変申し訳ございません。
詳細な回答いただき誠にありがとうございます。

当方で医療関係の統計処理の書籍を確認したところ、

1. 目的変数が正規分布している必要がある。
2. 目的変数、説明変数とも正規分布している必要がある。

という記載の本が非常に多く見られました。
医療系の人間は統計をあくまでツールとしてしか使用していないので、
学会でも上記の説が主流となっているようです。
正しい記述の本も、正規分布している必要がないものをわざわざ「正規分布している必要はない」と記載するわけではないので、上記の説が蔓延していると理解しました。

2値化してロジスティック回帰すると明らかに情報量の損失が発生しますので、
青木先生にご教授頂いたことを用いて、反論をしてみます。

ご丁寧にご教授頂き誠にありがとうございました。


計算関係にあるデータ間の相関検定の可否
  投稿者:コダマ 2019/08/25(Sun) 17:48 No. 22816
身長(A)と足の長さ(B)のデータが300人分あります。
胴長度(C)=1-B/Aで計算されるとき
BとCの相関を一次回帰分析で検定することは妥当でしょうか

背が高いほど胴長率が大きいのかを知りたいと思っているのですが、考えてみたら、ひとつひとつのデータ自体は一次計算式の代入変数と計算結果という関係にあるわけです。検定の可否をご教授いただけると幸いです。

Re: 計算関係にあるデータ間の相関検定の可否
  投稿者:青木繁伸 2019/08/25(Sun) 20:57 No. 22820
> ひとつひとつのデータ自体は一次計算式の代入変数と計算結果という関係

データ変換が一次変換であれば同じですが,この場合は一次変換ではないので,結果は異なるでしょう。
layout(1:2)
A = c(150, 160, 180, 175)
B = c(80, 75, 95, 82)
plot(B ~ A)
C = 1-B/A
plot(B, C)
layout(1)
predict(lm(B ~ A))
predict(lm(B ~ C))
の結果
> layout(1:2)
> A = c(150, 160, 180, 175)
> B = c(80, 75, 95, 82)
> plot(B ~ A)
> C = 1-B/A
> plot(B, C)
> layout(1)
> predict(lm(B ~ A))
1 2 3 4
75.71429 80.19780 89.16484 86.92308
> predict(lm(B ~ C))
1 2 3 4
87.51340 78.87029 86.76991 78.84639

Re: 計算関係にあるデータ間の相関検定の可否
  投稿者:コダマ 2019/08/26(Mon) 18:06 No. 22825
ありがとうございました。大いに参考にさせていただきます。


どう扱えば良いか困っています
  投稿者:りゅうへい 2019/08/25(Sun) 12:53 No. 22814
※答えを求めているわけでは断じてありません。ご理解下さい。
以下のデータが得られたとして 行う分析とそれにより得られる結論を述べよというものです。

会社A 資本金1000万 従業員100名
業績平均値:総務部70 営業部50 開発部50

会社B 資本金3億 従業員400名
業績平均値:総務60 営業30 開発60
いずれの平均値に対してもSD=5とする

分散分析を行うのだとは思うのですが,会社間では資本金と人数が異なるので
それをどう分析の中で考慮すればいいのかが分かりません。

ご教授のほどよろしくお願いします。

Re: どう扱えば良いか困っています
  投稿者:青木繁伸 2019/08/25(Sun) 17:25 No. 22815
業績平均値は各部の部員の業績の平均値ということでしょうかね?
資本金の違いは関係ないのでは?
全従業員数の違いも関係ないと思いますが,総務部,営業部,開発部 の人数がわからないので,分析不可能では?
つまり,会社全体の業績平均値を知ろうにも,各部の人数がわからなければ計算しようがない。

「これだけのデータではどのような分析もできない。問題自体が不適切という結論が得られる」が答えかな?

Re: どう扱えば良いか困っています
  投稿者:りゅうへい 2019/08/25(Sun) 18:36 No. 22818
ご回答ありがとうございます。
問題には 仮想の業績得点の平均値と書かれていました。
平均値の認識は ご指摘の通りで合っていると思います。
3部署合計人数=従業員数なのか 3部署以外の職員を含め 従業員数なのかが記載されていませんでした。 
3部署合計が全従業員数 という前提であれば 各部の人数が分からずとも
会社全体の平均値は Aは (70+50+50)/3=56.67 B (60+30+60)/3=50 
と 求めることは可能でしょうか?

資本金の違いは 結論が出た際の背景要因として触れる項目なのでしょうか?

ご教授のほどよろしくお願いします。

Re: どう扱えば良いか困っています
  投稿者:青木繁伸 2019/08/25(Sun) 20:40 No. 22819
> 3部署合計が全従業員数 という前提であれば 各部の人数が分からずとも会社全体の平均値は Aは (70+50+50)/3=56.67 B (60+30+60)/3=50 と 求めることは可能でしょうか?

不可能です。平均値の意味をわかっていれば,不可能であることは自明です。

たとえば,会社 A の合計人数が 100 人であっても,
総務部,営業部,開発部がそれぞれ,20, 30, 50 の場合は (70*20 + 50*30 + 50*50) / (20+30+50) = 54 だし,
22, 33, 45 の場合なら (70*22 + 50*33 + 50*45) / (22+33+45) = 54.4 だし,
1, 2, 47 の場合なら (70*1 + 50*2 + 50*97) / (1+2+97) = 50.2 だ。
こんな例を挙げる必要があるとは思わなかった(^_^;)

(70+50+50)/3=56.67 が成り立つのは, 総務部,営業部,開発部がそれぞれ,100/3, 100/3, 100/3 と等しいとき(あり得ないが)

(70*100/3 + 50*100/3 + 50*100/3)/(100/3 + 100/3 + 100/3) = 56.66666666

Re: どう扱えば良いか困っています
  投稿者:青木繁伸 2019/08/25(Sun) 21:00 No. 22821
> 資本金の違いは 結論が出た際の背景要因として触れる項目なのでしょうか?

たとえ他の条件が適切に示されていたとしても,資本金の違いがこのような違いをもたらしたかどうかについては,情報が「全く不十分」なので,何ともいえません。

この問題,インターネット上で提示されてるものですか?それなら,URL を示してください。

Re: どう扱えば良いか困っています
  投稿者:りゅうへい 2019/08/25(Sun) 21:15 No. 22823
すみません データ数が等しいケースだと勝手に解釈しておりました。
宿題で出された問題です。


どのような統計をつかうのでしょうか
  投稿者:こういち 2019/08/18(Sun) 03:28 No. 22802
青木先生
以下の疑問に関して教えてください

「疾患A型の患者」および「疾患B型の患者」における左室形態(4つタイプがありNormal, concentric remodeling(CR), concentric hypertrophy(CH), eccentric hypertrophy(EH))の分布を明らかにしました。結果の全体像は表1のようになりました。

ここで私が示したいことはA型ではCRが有意に多く、B型ではCHが有意に多い、
ということです。

そのために私がどのようにしたかというと
表2、表3のような2x2表をつくって(この作り方が正しいのか疑問ですが)
Fisher法で計算をしてそれぞれp=0.0427, p=0.0378という結果を得ました。
まず質問はこの統計処理が正しいのかを教えてください。

次にこれがもし正しいとしたら、下位比較をするという意味で、たとえば有意水準をBonferroniによって0.05/4C2=0.0083以下、とするのでしょうか?

別の考え方として4x2のクロス集計をおこなってP=0.1169を得たとして、
有意差なしでA型とB型の左室の型の分布に差がなし、で終了でしょうか?

以上、教えてください
よろしくお願いいたします


Re: どのような統計をつかうのでしょうか
  投稿者:青木繁伸 2019/08/19(Mon) 10:11 No. 22803
最初から CR,CH 型の比較だけを行うと決めていたのなら,(CR, CH)×(A,B) の 2×2 分割表の検定をすればよいのですが,今回はそうではなさそうなので,全ての組み合わせでの 2×2 分割表の検定をしてボンフェローニ法を適用する。(CR, CR 以外)×(A,B),(CH, CH以外)×(A,B)などとするのではないので注意。つまり,NormalとCR, NormalとCH,NormalとEH, CRとCH,CRとEH,CHとEH の 6 通りの比較を行う。
なお,「全体で有意である場合に下位検定として行える」という条件が必要かどうかは諸説ありますが,「この条件は不要」との見解が優勢かもしれない。

Fisher の正確検定でなければ,http://aoki2.si.gunma-u.ac.jp/lecture/Hiritu/Pmul-Tukey.html に書いてある方法を使えるでしょう(http://aoki2.si.gunma-u.ac.jp/lecture/Hiritu/Pmul-Ryan.html でもよいが,前者の法が通りがよいかな?)。

Re: どのような統計をつかうのでしょうか
  投稿者:こういち 2019/08/20(Tue) 00:56 No. 22804
青木先生

お返事をいただきありがとうございました。
responseが遅くなりお詫び申し上げます。
先生からいただいた説明を何度も何度も読み直してみましたが
私の説明が悪くて質問の意図が伝えられていないと思いました

私の疑問は
「たとえばA型という疾患においてCR, CH型の頻度を比較をする」ということではなく
「4つある形態の一つであるCRという形態を呈する患者が、A型では全体の56%, B型では全体の38%で、これらを比較するときに、どの様な統計を用いれば、CRという形態はA型がB型よりも有意に多いことを示せるか?」というA型とB型の比較のための統計処理がわからない、ということです。

そして同様に上記と独立して、
「4つある形態の一つであるCHという形態を呈する患者が、A型では全体の34%, B型では全体の52%で、これらを比較するときに、どの様な統計を用いれば、CHという形態はB型がA型よりも有意に多いことを示せるか?」というA型とB型の比較ののための統計処理がわからない、ということです。

私のいっている比較そのものがそもそも統計処理にそぐわず、有意に多い、ことは示せないのでしょうか?。
基本的なことが分かっていないとんちんかんな質問かもしれないと、
だんだん心配になってきましたが、どうか教えていただければ、本当にありがたいです。
よろしくお願いいたします。

Re: どのような統計をつかうのでしょうか
  投稿者:青木繁伸 2019/08/20(Tue) 08:28 No. 22805
No. 22802
> 私が示したいことはA型ではCRが有意に多く、B型ではCHが有意に多い、ということです。



No. 22804
> 「たとえばA型という疾患においてCR, CH型の頻度を比較をする」ということではなく...CRという形態はA型がB型よりも有意に多いことを示せるか

相反する説明がありますが。まあ,それは行方向のパーセントと列方向のパーセントのどちらを評価対象にするかということで,どちらが原因で,どちらが結果なのかで自ずと決まるものですが。でも,chisq.test でも,Fisher.test でも行・列の違いは無視されます(結果は同じ)。

> 「4つある形態の一つであるCRという形態を呈する患者が、A型では全体の56%, B型では全体の38%で、これらを比較するときに、どの様な統計を用いれば、CRという形態はA型がB型よりも有意に多いことを示せるか?」というA型とB型の比較のための統計処理

ならば,二項検定でよいでしょう。これだと,行ごとの検定になります。CR が A:66, B:19 で帰無仮説では母比率=0.5
> binom.test(66, 85)

Exact binomial test

data: 66 and 85
number of successes = 66, number of trials = 85, p-value = 3.041e-07
同様に CH では

> binom.test(40,66)

Exact binomial test

data: 40 and 66
number of successes = 40, number of trials = 66, p-value = 0.1089
ただし,この二つの検定は独立ではないです(独立性の検定という意味ではない)。
CR の 66 も 29 も Normal, CH, EH の度数による制約があります。

(Normal,CR,CH,EH)x(A,B) はそれらの制約のもと,Normal, CR, CH, H それぞれで,A, B の割合に違いがあるかを検定するものですよね。
検定ファミリーの一つに (CR,CH)x(A,B) があり,それは,CR, CH で A, B の割合に違いがあるかを検定する。
> fisher.test(matrix(c(66, 40, 19, 26), 2))

Fisher's Exact Test for Count Data

data: matrix(c(66, 40, 19, 26), 2)
p-value = 0.03111
CR では A,CH では B が有意に多い(A では CR,B では CH が有意に多い)といえるわけでしょう。

なお,また最初に戻りますが,「A型ではCRが有意に多く、B型ではCHが有意に多い」ということなら,A, B ごとに一様性の検定
> chisq.test(c(7, 66, 40, 5))

Chi-squared test for given probabilities

data: c(7, 66, 40, 5)
X-squared = 86.407, df = 3, p-value < 2.2e-16

> chisq.test(c(2, 19, 26, 3))

Chi-squared test for given probabilities

data: c(2, 19, 26, 3)
X-squared = 34, df = 3, p-value = 1.981e-07
でも,これはこの場合は Normal と EH がCR, CH より格段に少ないことが原因。
> binom.test(66, 106)

Exact binomial test

data: 66 and 106
number of successes = 66, number of trials = 106, p-value = 0.01478

> binom.test(19, 45)

Exact binomial test

data: 19 and 45
number of successes = 19, number of trials = 45, p-value = 0.3713
こんないろいろなことを考えなきゃいけないけど (CR,CH)x(A,B) でやればよいのでは?

Re: どのような統計をつかうのでしょうか
  投稿者:こういち 2019/08/21(Wed) 21:43 No. 22806
青木先生

丁寧に教えていただきまして、本当にありがとうございます。
教えていただいたことが恐らくは十分に理解されていないと思いますので
時間をかけて考えてみようと思います。

最後にひとつだけ教えてください、
先生が薦めてくださる「(CR, CH)x(A, B)をやる」、というのは
下のような2x2の表を作成してp=0.0311を求め、
(この表の作り方が正しいのかわかりません)
この場合「A型ではCRが多く、B型ではCHが多い」
(←この言い方が正しいかどうかにも自信がありません)ことが
0.0311>0.008(=0.05/6)であり有意とは言えない

という理解でよろしいでしょうか?
最後に教えていただければ本当にありがたいです


Re: どのような統計をつかうのでしょうか
  投稿者:青木繁伸 2019/08/22(Thu) 08:13 No. 22807
> A型ではCRが多く、B型ではCHが多い
> 0.0311>0.008(=0.05/6)であり有意とは言えない

でよいと思います。

Re: どのような統計をつかうのでしょうか
  投稿者:こういち 2019/08/22(Thu) 18:18 No. 22811
青木先生

お忙しい中いろいろ教えていただき
本当にありがとうございました
また、よろしくお願いいたします


3Dグラフの作成方法について
  投稿者:さくらとお城 2019/08/09(Fri) 14:19 No. 22791
青木先生

 いつもお世話になっております。

 複数のチャンネルで計測した経時データを,plot3Dなどで図示しようとしているが,なかなかうまく行きません。

 データとエクセルで作成したグラフは添付ファイルのようになります。実際は,もっとチャンネルが多いので,ここでは2つだけ表示しています。

 きれいに図示できる方法を教えてください。

Download:22791.zip 22791.zip


Re: 3Dグラフの作成方法について
  投稿者:青木繁伸 2019/08/09(Fri) 20:20 No. 22792
3次元グラフはあまりよい表示方法ではないでしょう。
添付する図を見ても,二次元グラフ(折れ線グラフ)のほうが二つのデータの相違・同等がよく分かるでしょう。
確かにたくさん重なるとわかりにくいでしょうが,3次元折れ線グラフが重なったら,なにがなんだかわからなくなります。


Re: 3Dグラフの作成方法について
  投稿者:青木繁伸 2019/08/09(Fri) 20:30 No. 22793
ggplot の例によくあるような,縦にずらっと折れ線グラフを並べるのはどうですか?
一つのグラフに,傾向の似通った複数の折れ線グラフを描けばなおよいかも。


Re: 3Dグラフの作成方法について
  投稿者:さくらとお城 2019/08/10(Sat) 10:56 No. 22794
青木先生

 ご返事ありがとうございます。添付のようなグラフを考えていました。ちなみに,これが別ソフトで作ったもので,細かい設定などが大変で,Rでもできるかなと思って質問しました。

 先生のggplotの例も良さそうで,試してみます。


Re: 3Dグラフの作成方法について
  投稿者:さくらとお城 2019/08/22(Thu) 10:22 No. 22808
青木先生

 その後,色々と調べたけど,恥ずかしながら,なかなか先生のような図を描けません。ggplotの例を書いたコードを見せて頂いて宜しいでしょうか?

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

Re: 3Dグラフの作成方法について
  投稿者:青木繁伸 2019/08/22(Thu) 12:50 No. 22809
私自身,ggplot は使わないので,細かな設定は知りません。
以下のようなプログラム例を示しておきますので,追加,修正の上ご利用ください。
ポイントは通常のデータフレームを時間とチャンネルとデータ値の3列から成るデータフレームに変換して,それを用いて折れ線グラフを描くということです。
set.seed(123)

m = 10
n = 50

# n x m 行列:列ベクトルがチャンネルごとの時系列データ

simdata0 = matrix(0, n, m)
for (i in 1:m) {
simdata0[, i] = cumsum(rnorm(n))
}

# チャンネルごとのデータを,
# time:時間(1:n),value:データ値,チャンネル名の,n x 3 のデータフレームとして rbind で縦に繋いだデータフレームに変換

simdata = NULL
for (i in 1:m) {
simdata = rbind(simdata, data.frame(time = 1:n, value = simdata0[, i],
channel = paste0("x", i)))
}

g = ggplot(simdata, aes(x = time, y = value))
g = g + geom_line(aes(colour = channel), size = 0.5)
g = g + geom_point(aes(colour = channel), size = 1, alpha = 0.5)
g = g + facet_wrap(~channel, ncol = 2) # ncol で何列のグラフにするか
g = g + theme_classic() # または g = g + theme_bw() でも
g
以下のようなグラフになります(部分)


Re: 3Dグラフの作成方法について
  投稿者:さくらとお城 2019/08/22(Thu) 12:53 No. 22810
青木先生

 早速のご返事ありがとうございます。先生のコードを勉強させていただきます。


統計分析手法の選択につきまして
  投稿者:こだま 2019/08/14(Wed) 13:13 No. 22795
青木先生

はじめまして。

自身の周りに統計手法を専門的に理解している人間があまりおらず書き込みさせて
いただきます。

私は現在、大学院生(博士後期課程)に在学中で現在GISを用いた緑地の分布に
ついての分析結果学を会論文投稿に向け作業をしているところです。

その中で、現在の緑地の分布要因は何に規定されているのか?を統計手法による
分析を用いて明らかにしたいと考えております。

その際にどの分析手法が良いのかに行き詰っているところです(判別分析なのか、
数量化2類なのかと悩んでいます。)。

データはメッシュで集計しており,メッシュ毎に以下のデータを格納しています。

目的変数:緑地の分布タイプ(別途分析により7タイプに分類済)
説明変数:用途地域(13タイプ),区画整理事業の有無,傾斜角度,標高,駅の有無等

目的変数が質的データ、説明変数が質的データと量的データが混在しているため、単純
に判別分析や数量化2類を使えない?のかと思っています(理解不足なのかもしれませんんが)。

お忙しいところお手数をおかけしますが、どうぞご回答いただけますと幸いです。

Re: 統計分析手法の選択につきまして
  投稿者:青木繁伸 2019/08/14(Wed) 13:20 No. 22796
> 目的変数が質的データ、説明変数が質的データと量的データが混在しているため、単純
に判別分析や数量化2類を使えない?のかと思っています

の点についてのみ。
質的データはダミー変数として取り扱うことで,量的データと一緒にして判別分析に使用できます。(ちなみに,数量化II類はダミー変数を使う判別分析と等価です)

Re: 統計分析手法の選択につきまして
  投稿者:こだま 2019/08/14(Wed) 13:58 No. 22797
青木先生

早速のお返事ありがとうございます。

判別分析を採用してみたいと思います。

標高等の量的なデータを区分して数量化2類かな?などと考えて
いて、その区分の方法が難しいなと思っていたところでした。

今後、また相談させていただくことがあるかもしれませんが、
どうぞよろしくお願いいたします。

Re: 統計分析手法の選択につきまして
  投稿者:こだま 2019/08/14(Wed) 17:12 No. 22798
青木先生

度々申し訳ございません。

色々な文献(多変量解析の本やインターネット、判別分析を行っている研究論文)を
調べているのですが、多くが目的変数が2択(例えば、公園がある/ない のような者)
の事例でした。

今回、分析しようとしているように目的変数が7つのタイプがあっても分析手法として
採用しても問題はないのでしょうか。

可能であればで結構ですので、ご教授いただけますと幸いです。

Re: 統計分析手法の選択につきまして
  投稿者:青木繁伸 2019/08/14(Wed) 21:32 No. 22799
ロジスティック判別も線形判別も,3群以上であっても判別はできます。
どこでもよいですが,判別分析のページを読んでみてください。
一番簡単な2群判別の話しか書いていないのは,くそページです。

7群判別ともなると理論的には可能ですが,実際に分析して得られる解の解釈は難しいかも。(貴方の用意した説明変数で十分に判別できるかどうかということですが)

Re: 統計分析手法の選択につきまして
  投稿者:こだま 2019/08/14(Wed) 23:06 No. 22800
青木先生

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

論文提出まで時間があるので判別分析含め、その他の多変量解析について
色々な参考書やページを読んで勉強いたします。

説明変数の選択が難しいところですが、こちらは既往文献等を参考にしつつ
選択していきたいと思います。

実際に分析を進めてみて何かありましたら相談させていただきたいと思います。

夜分遅くまでお手数をおかけしました。


【R】リストに格納された複数データフレームのEXCEL保存
  投稿者:明石 2019/08/05(Mon) 16:12 No. 22787
青木先生 様;

お忙しいところを失礼いたします、明石と申します。
毎々、ご丁寧なご教示をいただき、誠にありがとうございます。
改めて、御礼を申し上げます。

青木先生にご教示いただきたいことが出てきました。
何卒どうぞよろしくお願いいたします。

----------------------------------------

リストに格納された複数データフレームを、
拡張パッケージ openxlsx を用いて、
EXCELの複数のシートに格納したxlsxファイルで保存することを考えます。

以下、irisデータでお示しをします。

data(iris)

dfs <- split(iris, iris$Species)

library(openxlsx)  # 拡張パッケージ

write.xlsx(dfs, file="iris_Species.xlsx")

setosa、versicolor、virginica の3つのシートからなる
"iris_Species.xlsx"ファイルは保存されますが、以下のメッセージがでます。

「Note: zip::zip() is deprecated, please use zip::zipr() instead」

私が書いたコマンドには、zip()関数はありませんので、
このようなメッセージが表示された場合、どのように対処すればよいのでしょうか。

ご教示をいただけましたら、大変に助かります。
お手数をおかけいたします。

Re: 【R】リストに格納された複数データフレームのEXCEL保存
  投稿者:青木繁伸 2019/08/05(Mon) 19:40 No. 22788
write.xlsx から呼ばれる saveWorkbook の中で zip を呼ぶようですが(明示的には見えない),私にはそれ以上追求できませんでした。
deprecated というのは,将来的に廃止されるということで,無効になっているわけではないです。ということで,しばらくはそのまま使えば良いと思います。そのうち,使えなくなったらまずは openxlsx パッケージのバージョンを確認して,新しくなっていたらインストールし直すということになるでしょう。それもできないという事態になったら,パッケージ開発者に連絡するしかなくなるかな。

Re: 【R】リストに格納された複数データフレームのEXCEL保存
  投稿者:青木繁伸 2019/08/05(Mon) 19:43 No. 22789
Excel シートにする必要性がどの程度あるのか,他の人と共有するためには Excel シートでなければならないとか,いきなり Excel を書き込むのではなく,CSV ファイルとして書き出すのではだめなのか,とか,私なら Excel を使わない方向で行きたいですね。

Re: 【R】リストに格納された複数データフレームのEXCEL保存
  投稿者:明石 2019/08/06(Tue) 08:05 No. 22790
青木先生 様;

お忙しいところを失礼いたします、明石と申します。
毎々、ご丁寧なご教示をいただき、誠にありがとうございます。

今回も有難いご教示をいただき、大変に助かりました。

青木先生がお書きになられたことを、戒めたいと思います

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


【R】リストに格納された複数データフレームのファイル保存(ループを使わない)
  投稿者:明石 2019/07/29(Mon) 17:25 No. 22784
青木先生 様;

お忙しいところを失礼いたします、明石と申します。
毎々、ご丁寧なご教示をいただき、誠にありがとうございます。
改めて、御礼を申し上げます。

青木先生にご教示いただきたいことが出てきました。
何卒どうぞよろしくお願いいたします。

--------------------------------------------------

リストに格納された複数データフレームのファイル保存を、
ループを使わない書き方について、ご教示いただければ大変に助かります。

以下、例示いたします。

データフレーム dat を、カラム名"支店名"でファイル分割したとします。

dfs <- split(dat, dat$支店名)

支店名ごとに分割されたデータフレーム群を、csvファイル保存します。

n <- length(dfs)
names <- names(table(dat$支店名))

for (i in 1:n) {
 fn <- paste(names[i],".csv", sep="")
 write.csv(dfs[[i]],fn, row.names=F)
}

ループを回さないで、apply関数群などを用いて、シンプルに書けるようならば、
ご教示をいただきたいと思います。

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

Re: 【R】リストに格納された複数データフレームのファイル保存(ループを使わない)
  投稿者:青木繁伸 2019/07/29(Mon) 22:16 No. 22785
iris データセットでやってみましょう

Species で 3 分割

iris2 <- split(iris, iris$Species)

書き出す。junk への代入は,ゴミを表示させないためだけ。

junk <- sapply(names(iris2), function(n) write.csv(iris2[[n]], paste(n, ".csv", sep="")))

あまりわかりやすいものではないですね。
貴方が示した for ループを使うほうがずっと分かりやすい。

Re: 【R】リストに格納された複数データフレームのファイル保存(ループを使わない)
  投稿者:明石 2019/07/30(Tue) 06:13 No. 22786
青木先生 様;

お忙しいところを失礼いたします、明石と申します。
毎々、ご丁寧なご教示をいただき、誠にありがとうございます。

今回も、大変に良い勉強をさせていただきました。

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


因子分析の因子の差について
  投稿者:コロン 2019/07/16(Tue) 13:38 No. 22778
いつもお世話になっております。

早速ではございますが,「4つの因子」が出てきたという例で説明させて頂きます。

数学力が高いもの,中程度のもの,低いものという3つのグループの間で,因子において差があるかどうかを見る場合なのですが...

(1)因子毎に,一要因分散分析(被験者間)を行い,差があれば多重比較で見ていく
(2)4つの因子をすべて用いて,二要因分散分析(被験者「間」×被験者「内」)を行い,必要に応じて下位検定を行う

どちらが正しいのでしょうか?古い論文を見ていると(1)が圧倒的に多いのですが...。

よろしくお願いいたします。

Re: 因子分析の因子の差について
  投稿者:青木繁伸 2019/07/17(Wed) 09:06 No. 22779
理論的には (2) が正しい(妥当)と思います

なぜ (1) が圧倒的に多いかは,いくつか理由があるでしょう

サンプルサイズが小さく,(2) に耐えない
モデルが複雑になり,解釈が難しくなる
コンピュータプログラムがなく(利用しにくく),計算できない
「だって,みんな (1) をやっているんだもの」と長いものに巻かれる
そもそも (2) という選択肢を知らない

Re: 因子分析の因子の差について
  投稿者:コロン 2019/07/17(Wed) 09:36 No. 22780
青木先生

お忙しい中,いつもありがとうございます。

やはり,そうですよね....。かなり検索をしてみたのですが,すべて(1)で処理されていたのです。(2)が適切ではないかなとずっと考えていましたので,すっきり致しました。

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

Re: 因子分析の因子の差について
  投稿者:コロン 2019/07/20(Sat) 17:19 No. 22781
青木先生

もう一つ重要なことを聞き忘れていました。

各因子の項目の数が異なるのですが、尺度得点あるいは因子得点を利用して2要因分散分析を行ってもよろしいでしょうか。

よろしくお願いいたします。

Re: 因子分析の因子の差について
  投稿者:青木繁伸 2019/07/24(Wed) 11:04 No. 22782
項目の数が異なるというのが二元配置分散分析の障害になるのですか?

Re: 因子分析の因子の差について
  投稿者:コロン 2019/07/25(Thu) 16:34 No. 22783
青木先生

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

【要因1】
因子1 7項目
因子2 5項目
因子3 4項目

【要因2】
学力上位群
学力下位群

項目数が異なるため,各因子の平均値をそのまま解釈して良いのかなと悩んでおりました。ひょっとしたら,過去の研究,項目数が異なるために因子毎に比較していたのかなと考えたり...。

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


固辞,数字の置換
  投稿者:コロン 2019/07/08(Mon) 13:44 No. 22773
いつもお世話になっております。

基本的な質問で申し訳ございませんが,以下のデータがあるとします。

class score
1 10
1 11
2 13
2 14
3 12
3 15
4 17
4 14
5 13
5 10

classの 3 を 5 に修正して新しいデータフレーム形式で表示させるにはどのようにすればよろしいのでしょうか?

Re: 固辞,数字の置換
  投稿者:青木繁伸 2019/07/08(Mon) 15:37 No. 22774
データフレームを a とすれば,
a$class[a$class==3] <- 5
で置換できます。

Re: 固辞,数字の置換
  投稿者:コロン 2019/07/08(Mon) 16:40 No. 22775
青木先生

お忙しい中申し訳ございません。

ネットで検索をしていたら,置換は gsub というサイトが多かったため,gsubでやっていたのですが,何度やっても,データフレーム形式にならず,かつ,すべての数字が5に変換されたので,わからなくなり,投稿した次第です。

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

Re: 固辞,数字の置換
  投稿者:青木繁伸 2019/07/08(Mon) 18:07 No. 22776
gsub でもできないことはないですが,gsub は本来は文字列の置換ですから
a$class <- as.numeric(gsub(3, 5, a$class))
as.numeric してやるなんて余計なことをしないといけません。

Re: 固辞,数字の置換
  投稿者:コロン 2019/07/08(Mon) 18:24 No. 22777
青木先生

as.numericなんですね!gsubでもできました。ただ先生から最初に提示いただきましたものがわかりやすいです。

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


hclust()の距離行列とメソッドの指定
  投稿者:波音 2019/06/28(Fri) 15:03 No. 22764
Rのhclust()について教えてください。

引数にward.Dとward.D2が指定できるようになっていますが、例えば「d」という距離行列(ユークリッド距離)が存在したと時にいくつかの方法がネット上や書籍で示されています。

・距離行列をそのまま指定する
hclust(d, method="ward.D")

・距離行列の自乗したものを指定する
hclust(d^2, method="ward.D")

・川端ほか「Rにおる多変量解析入門」オーム社(2019)P374記載
hclust((1/2)*d^2, method="ward.D")

結局のところ「解釈ができれば何でもいい」というくらいの問題なのかもしれませんが、アルゴリズム的(?)にはどういう指定が適切なのでしょうか。

川端ほか「Rにおる多変量解析入門」オーム社(2019)P313にward.Dとward.D2の違いが説明されていますが、どういう距離行列とメソッドの組み合わせ指定をすればよいのでしょうか、、、

Re: hclust()の距離行列とメソッドの指定
  投稿者:青木繁伸 2019/06/29(Sat) 12:05 No. 22765
川端さんが具体的にどのように説明されているのか分かりませんが,私の理解は
http://aoki2.si.gunma-u.ac.jp/lecture/misc/clustan.html
に示してあります。これは奥野忠一らの多変量解析,多変量解析入門にあったものだと思います。まだパソコンなんてものがなく,コンピュータはFORTRAN,COBOLでプログラムするという時代のものです。もともとは,CLUSTAN というプログラムのアルゴリズムの説明にあるものです。実際,私もこれに基づいてクラスター分析のプログラムを書いて使っていました。
当初の R には今の ward.D が ward として実装されていました。推測するに,「バグ」だったのだと思います。
私の書いた本でも,「ウォード法と重心法を採用する場合にはユークリッド平方距離を使うべきである。ユークリッド平方距離を選択できるオプションがないので,後述の hclust の第 1 引数としてユークリッド距離を二乗したものを指定することを忘れないようにしよう。」と書きました。
その後指摘があったかなんかで,本来あるべき ward.D2 が実装されましたが,過去のしがらみを捨てきれなくて(バグだったとも言えなくて)両方を残したのだと邪推しています。

Re: hclust()の距離行列とメソッドの指定
  投稿者:波音 2019/07/01(Mon) 10:41 No. 22769
ご回答ありがとうございます。

たしかに青木先生の本にも「hclust の第 1 引数としてユークリッド距離を二乗したものを指定」と書かれており、それ以前にもこのように教わっていたので、私の中ではある意味ルール的にd^2としていましたが、、、

やはり(使う人の立場によっては程度問題ではあると思いますが)アルゴリズムを理解しておくことは重要だと感じさせられました。

少し話はそれてしまいますが、誰でも使えるものというのは「考えなくても(理解していなくても)使えるもの」ともいえるので、ただでさえ色々と誤解の多い統計話題は今後さらに注意しないといけないと思いました。

※パッケージ文化に恵まれた我々は特に注意せねばならないのかもしれません。


ダミー変数を従属変数にした場合の切片
  投稿者:山本 2019/06/29(Sat) 23:47 No. 22766
はじめて書き込みます。山本です。
講義の関係で重回帰分析を行なっているのですが、ダミー変数を従属変数とした際に切片が負という結果がでました。
ダミー変数は0か1かですので正になるはずなのですが…なぜでしょうか
この結果はどう解釈すればよいのでしょうか
それとも説明変数を変更する必要があるのでしょうか
お願いいたします

Re: ダミー変数を従属変数にした場合の切片
  投稿者:青木繁伸 2019/06/30(Sun) 07:19 No. 22767
従属変数が 0/1 の二値データであることと,切片が負の値を取ることは何の矛盾もありません。
というか,独立変数のとる値により切片は負の値になることも,[0, 1] の値になることも,1より大きな値になることもあります。
たとえば,
d <- data.frame(y = c(1, 0, 1, 1, 0, 0, 0, 1, 0, 0),
x1 = c(58, 37, 38, 50, 36, 59, 49, 66, 47, 58),
x2 = c(54, 53, 50, 46, 46, 65, 37, 69, 37, 44))
の場合は,切片は負の値になります。
> lm(y~x1+x2, data=d)
(Intercept) x1 x2
-0.646540 0.005094 0.015825
ここで,x1 から 250 を引いた変数 x1.250 を作り,それを用いて予測式を計算すると
> d$x1.250 <- d$x1 - 250
> lm(y~x1.250+x2, data=d)
(Intercept) x1.250 x2
0.627073 0.005094 0.015825
ここで,x1 から 435 を引いた変数 x1.435 を作り,それを用いて予測式を計算すると
> d$x1.435 <- d$x1 - 435
> lm(y~x1.435+x2, data=d)
(Intercept) x1.435 x2
1.569546 0.005094 0.015825
それぞれの偏回帰係数は全く同じになっていることに注意。
つまり,独立変数の平行移動をした場合,その影響を打ち消すために定数項が調整されるということ(ちなみに独立変数を定数倍したときは偏回帰係数が定数分の1になる)。

で,今回の問題の本質は,なぜ従属変数が 0/1 の二値データなのに予測値が0より小さくなったり,1より大きくなったりしてしまうのかということ。

それは直線回帰をするからやむを得ないということです。
図を見れば明らかです。

二値データはロジスティック回帰分析を行うのが定石です。
図の赤で示したのがロジスティック曲線です。
> glm(y~x1+x2, data=d, family=binomial(link="logit"))

Call: glm(formula = y ~ x1 + x2, family = binomial(link = "logit"),
data = d)

Coefficients:
(Intercept) x1 x2
-5.08999 0.02083 0.07214

Degrees of Freedom: 9 Total (i.e. Null); 7 Residual
Null Deviance: 13.46
Residual Deviance: 11.93 AIC: 17.93


Re: ダミー変数を従属変数にした場合の切片
  投稿者:山本 2019/06/30(Sun) 12:02 No. 22768
素早く回答していただきありがとうございます。
そういった分析もあったのですね…
よく分かりました


介入研究の因子分析
  投稿者: 2019/06/13(Thu) 00:34 No. 22746
青木先生

いつもお世話になっております。
また、教えていただきたいことがあり、ご連絡させていただきました。

現在、先行研究で使われていた尺度(20項目ほど)を用いて調査を行いました。
その先行研究では因子分析を行っていませんでしたが、今回の調査は介入研究なので、尺度を因子分析を行った方が傾向を示しやすいと思うのですが、介入研究の場合の因子分析は、どのように行えば良いのでしょうか?
ちなみに、介入前と介入後の得点で、それぞれ因子分析しましたが、同じように分かれませんでした。

ご助言いただけますでしょうか。
よろしくお願いいたします。

Re: 介入研究の因子分析
  投稿者: 2019/06/17(Mon) 15:24 No. 22758
青木先生

いつもお世話になっております。

先日、6月13日にご質問させていただいた者です。
どうも、タイミング悪く入ってしまったようで、お返事を、まだ、いただけておりません。
ぜひとも、ご助言いただけないでしょうか?

よろしくお願いいたします。

Re: 介入研究の因子分析
  投稿者:青木繁伸 2019/06/17(Mon) 19:45 No. 22759
上の方にあるものばかり見ていたようで,気づきませんでした。

さて,先行研究で使われていた尺度は因子分析で作成されたものではないと言うことでしょうか。

とすると,まずいですね。あなたの調査データで一因子性がないという結果になっても,先行研究を否定するということになってしまいますし,対象が違うから違う因子が得られたのだといわれるかも知れませんし。

介入研究だからといって,因子分析のやり方に特徴や注意点があるとは思えませんが。

介入前後で因子が異なることは想像するに難くないと思いますよ。

実際の因子分析の結果だけでなく,基礎的な統計解析を見ずにああだこうだいうのは難しいです。もう少しはっきりいえば,あなたの持っているデータを自由にいろいろ分析してみてその結果を総合して何かをいうということでなければ(貴方の分析結果を見ただけではなく),なんもいえない。ということです。そういうことをするのは,共同研究者以上の立場ですよね。査読者はある意味そのような立場にあるのでしょうが,実際に「データを見せろ,私が解析してみて貴方の言うことが正しいかどうか示してやる」とは言えない訳です。

データ解析というものは,そういうものです。

Re: 介入研究の因子分析
  投稿者:青木繁伸 2019/06/17(Mon) 21:33 No. 22760
「十分な分析を行い,以下のような結果を得ました」というとき,どのように十分な分析だったかを「詳細に」示し,それを査読者が十分だったと認めるかどうかということです。

蛇足ではありますが,査読者が十分だったと認め論文掲載となっても,読者(第三者)の十分だったかどうかの判断とは必ずしも一致しない(投稿論文への異議申し立てはあり得る)。

その昔,ある学会誌に「奇病」の投稿があったが,老臨床医が「それは,脚気だろ」と指摘してチャンチャン!となったことがある(そんなに大昔のことではない)

厳しく言えばそういうこと。

Re: 介入研究の因子分析
  投稿者: 2019/06/18(Tue) 00:41 No. 22761
青木先生

ご助言ありがとうございます。

そもそも、因子分析されていない尺度を使う事自体、無理があるということですね。
もう1つ、気がかりなことがあります。Nの数が30に満たないので、それも、因子分析するには、少なすぎるかなあと思っておりました。

このまま、因子分析を行うことは、諦めた方が良いのかとも思いますが、以下のように、進めるのは、間違いでしょうか?

介入前の因子分析の結果では、4因子に分かれ、α係数も1つだけ0.7に足りませんでしたが、3つは0.8以上あり、内的整合性は高いと判断できます。
介入後も、介入前に抽出された因子に合わせて、α係数を出して充分な値であれば、その後に介入前と介入後で対応のあるt検定をする。

・・・というのは、分析としてはどうでしょうか?
ご意見を伺えると助かります。

Re: 介入研究の因子分析
  投稿者:青木繁伸 2019/06/18(Tue) 08:38 No. 22762
> Nの数が30に満たないので、それも、因子分析するには、少なすぎるかなあと思っておりました。

そこは最初に書いておいてくれなくては。
圧倒的に少なすぎます。後で述べる対応のある t 検定であっても,十分とは言いがたい。

> その後に介入前と介入後で対応のあるt検定をする。

何の対応のある t 検定ですか?先行研究の尺度の合計値ですか?それとも,4因子それぞれに含まれる項目の合計点ですか?

そもそも,因子分析の結果の4因子はそれぞれ何項目が含まれるのか,それぞれの因子の意味づけ(解釈)はできるのか。

もともと20項目あるなら,項目別に検定を行ったらどうなるのか?

Re: 介入研究の因子分析
  投稿者: 2019/06/18(Tue) 12:13 No. 22763
青木先生

おはようございます。
ご回答いただいて、ありがとうございます。

項目別の検定は、既に行っていて、効果は確認済みです。
しかし、ご指摘いただいたように、Nの数が少ないので、今回は因子分析は諦めようと思います。
しかし、今後、Nを増やして再度、分析を行いたいと思っているので、このまま、続けさせていただきたいと思います。

対応のあるt検定については、抽出した各因子の下位尺度得点として、介入前、介入後の各平均値を算出し、比較をしようと思っています。

しかしながら、今後、因子分析を行うときには、各因子内の項目の数や、各因子の意味づけ(解釈)に気を付けていきたいと思います。
そして、因子分析を用いた方が良いのか、項目別に検定だけの方が良いのか判断しようと思います。

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


【R】do.call()関数を用いたリスト(複数のデータフレーム)のmerge
  投稿者:明石 2019/06/16(Sun) 09:50 No. 22754
青木先生 様;

お忙しいところを失礼いたします、明石と申します。
毎々、ご丁寧なご教示をいただき、誠にありがとうございます。
改めて、御礼を申し上げます。

青木先生にご教示いただきたいことが出てきました。
何卒どうぞよろしくお願いいたします。

--------------------------------------------------

リスト(複数のデータフレーム)のmergeを、forループで順次mergeするのではなく、
do.call()関数を用いて、もし書けるようならば、ご教示をお願いいたします。

do.call(rbind, list) にならい、
do.call(merge, list) のように、もし書けたとしたら、どのように書けるのか興味があります。

# tempフォルダ下にある、csvファイルのファイル名を取得

path <- setwd("C:/temp")
names <- list.files(path, pattern="csv")
nn <- length(names)

d <- list(nn)

loop <- 1: nn

for (i in loop) {
d[[i]] <- read.table(names[i], header=TRUE, sep=",", stringsAsFactors=F)
}

# リスト d に格納されている nn個のデータフレームは、いずれも同じ、mergeの属性 "支店名"をもっています。
# 要素の添字に沿って、完全外部結合を順次、繰り返したいと思います。

for (i in loop) {
if( i == 1) {
dat <- d[[1]]
} else {
dat <- merge(dat, d[[i]], by.x="支店名", by.y="支店名" , all.x = TRUE, all.y = TRUE ,no.dups = F)
}
}

もし、do.call(merge, d) のように書けるようならば、ご教示をいただきたいと思います。

お手数をおかけいたします。
何卒どうぞよろしくお願いいたします。

Re: 【R】do.call()関数を用いたリスト(複数のデータフレーム)のmerge
  投稿者:青木繁伸 2019/06/16(Sun) 12:34 No. 22756
別に do.call でなくても,簡単に書きたいということなんでしょうか?

Reduce を使うのはいかが?
set.seed(11111)
d <- vector("list", 3)
d[[1]] <- data.frame(a=1:5, b=letters[1:5])
d[[2]] <- data.frame(a=c(1,3,5,6,8), c=LETTERS[1:5])
d[[3]] <- data.frame(a=c(2,4,5,7,9), d=rnorm(5))

merge2 <- function(x, y) merge(x, y, all.x = TRUE, all.y = TRUE, no.dups = FALSE)
Reduce(merge2, d)
一行で書くと
Reduce(function(x, y) merge(x, y, all.x = TRUE, all.y = TRUE, no.dups = FALSE), d)
なお,for を使う場合も,以下のようにすれば十分簡潔に書けるのでは?
x <- d[[1]] # 結果の初期化は,最初のデータフレームを代入
for (i in 2:length(d)) { # 2 番目以降のデータフレームをマージ
x <- merge(x, d[[i]], all.x = TRUE, all.y = TRUE, no.dups = FALSE)
}
x
どのように書いても,以下のような結果になると思いますが。
  a    b    c            d
1 1 a A NA
2 2 b <NA> 0.003630422
3 3 c B NA
4 4 d <NA> 0.798122717
5 5 e C 0.893397342
6 6 <NA> D NA
7 7 <NA> <NA> -1.195893897
8 8 <NA> E NA
9 9 <NA> <NA> -1.269878611


Re: 【R】do.call()関数を用いたリスト(複数のデータフレーム)のmerge
  投稿者:明石 2019/06/16(Sun) 16:41 No. 22757
青木先生 様;

お忙しいところを失礼いたします、明石と申します。

先生がお書きになられたように、簡潔に書きたいというのが意図です。
そこで、do.call()関数に着目しました。

幾つかの方法をお示しくださいましたので、よく理解できました。
reduse()関数を初めて知りましたので、これを機会に勉強させていただきます。

今回も助けていただき、大変に助かりました。
誠にありがとうございました。


低値・高値の基準について
  投稿者:H 2019/06/14(Fri) 11:49 No. 22750
青木先生
はじめまして。

ホルモン剤内服後のホルモン値が妊娠率に影響するかについて調べております。
ホルモン値低値群・高値群に分け、その間のコントロール群と妊娠率を比較したいのですが、どのようにして低値・高値の基準を設定すればいいでしょうか。
他の研究でもホルモン剤内服後のホルモン値が妊娠率に影響するかについての報告がありますが低値と高値に分け検討しているものがなく検討したいと思っています。

他の研究と内服薬の種類が異なるためか他の研究とホルモンの平均値が異なっていたため、第1四分位数以下、第3四分位数以上で低値・高値を分けコントロール群と比較し有意差がでたのですが恣意的と言われました。

多忙のところ申し訳ありませんが、適切な基準の設定があればご教授お願い致します。

Re: 低値・高値の基準について
  投稿者:青木繁伸 2019/06/14(Fri) 21:04 No. 22751
> 他の研究でもホルモン剤内服後のホルモン値が妊娠率に影響するかについての報告がありますが

そのような報告ではどのような分析手法が使われていますか?その分析手法を踏襲することに何らかの困難があるのでしょうか?

どのような分析手法が使われているのか知らない状態ではありますが,私がおすすめできる分析手法は「ホルモン値を独立変数,妊娠の有無を従属変数として,ロジスティック回帰分析」をすれば良いのではないですか?
独立変数としてはホルモン値の他の要因も考慮できますし。

データはなるべくもとのまま使うのが原則でしょう。何らかの加工(第1四分位数以下、第3四分位数以上で低値・高値を分けコントロール群と比較し)などを行うと,必ず恣意的ではないか?といわれるのは避けがたいことです。では,といって,別の分類法も使って,別の結果が出たら,どっちが正しいのかといわれるし,いい結果が出たものだけを発表すると,結局は「恣意的だろう」といわれる。

Re: 低値・高値の基準について
  投稿者:H 2019/06/15(Sat) 01:33 No. 22752
返信ありがとうございます。
他の報告ではホルモン値が100毎に測定してあるものや平均値の低値と高値で評価してあるものがあります。ただそれらでは明確な根拠がないため恣意的になるかと思われます。

投降後にロジスティック回帰を行いましたが有意差を認めませんでした。
低値群で妊娠率が低いため、統計的に評価したいのですが恣意的になってしまいますか。
宜しくお願い致します

Re: 低値・高値の基準について
  投稿者:青木繁伸 2019/06/15(Sat) 13:15 No. 22753
> 投稿後にロジスティック回帰を行いましたが有意差を認めませんでした
ということなら,仕方ないでしょう。

有意ではないが医学的には意味があるようなら,奥の手で,もう少しデータを増やすとか。

Re: 低値・高値の基準について
  投稿者:H 2019/06/16(Sun) 11:50 No. 22755
返信ありがとうございます。
有意差なしでデータを出そうと思います。
ありがとうございました。


サンプルについて
  投稿者:林檎 2019/06/12(Wed) 17:25 No. 22744 HomePage
青木先生

はじめまして。
サンプルについてわからないことがあり、質問をさせていただきました。

私は、10施設の企業を対象に質問紙調査による横断研究を実施し、解析をしています。
500名の有効回答があり、施設ごとと従属変数にて一元配置分散分析を行ったところ、1施設のみ合計得点平均が有意に高いことが分かり、指導教授より、その施設は特殊であるため抜いたほうがいいのではないかとご指摘を受けました。ちなみにこの施設の対象は20名でした。
世の中の論文をみると、このような有意差を見ている論文はさまざまであり、今回先生のご意見をお聞かせいただけましたらと思っております。

宜しくお願い致します。

Re: サンプルについて
  投稿者:青木繁伸 2019/06/12(Wed) 22:03 No. 22745
> その施設は特殊であるため抜いたほうがいいのではないかとご指摘

そもそも,調査の目的は何だったのですか?というか,調査をして,何が分かると思っていたのですか?

「どの施設も同じようである」ということを知りたかったのなら,「1施設のみ合計得点平均が有意に高い」というのは,ビックリするような知見ではないですか。その施設は何らかの問題があると。

対象者が20名ということが特殊なのか,それ以外の条件が特殊なのか。いずれにしても,サンプルサイズが20であっても,それでも有意な差であると言うことですよ。

「その施設は特殊であるため抜いたほうがいい」ということなら,最初からその施設は対象としない方が良かったでしょう。不都合な結果が出たからなかったことにしようというのは,はなはだ科学的ではないです。

その施設を除いて,「どの施設も同じようでした」という結果を出して,何の意味があるのでしょうか?

指導教授とよく話し合うことをお勧めします。

Re: サンプルについて
  投稿者:林檎 2019/06/13(Thu) 08:02 No. 22747 HomePage
青木先生

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

目的としては、重回帰分析を行い、企業の社員のストレスの要因を明らかにすることです。
重回帰分析に入れるための変数を探すこと、属性を明らかにするために一元配置分散分析を実施したら、有意差があったという結果でした。
調査をお願いするときには、似たような企業だと思い調査を行ったのですが他の施設よりストレスが有意に高かったことが分かりました。

指導教授と相談をしてみます。
ありがとうございました。

Re: サンプルについて
  投稿者:青木繁伸 2019/06/13(Thu) 10:32 No. 22748
> 他の施設よりストレスが有意に高かったことが分かり

それこそ,その企業が他の企業とどこが違うか(ストレスの原因)を知るために使える重要なサンプルではないですか。

Re: サンプルについて
  投稿者:林檎 2019/06/13(Thu) 14:00 No. 22749 HomePage
青木先生

ご返答いただきありがとうございます。

この企業は大切なサンプルであるということを改めて実感致しましたので、先生からのお返事をもとに解析を進めていきたいと思います。

とても丁寧なご対応をいただき、ありがとうございました。


クラスタリングされた個体の割合の比較とライアンの方法について
  投稿者:PV 2019/06/09(Sun) 19:11 No. 22739
青木先生、

はじめまして、
ご多忙のところ失礼いたします。

クラスタリングされた個体を群間で比較し、割合が有意に異なるかどうか、またどの群間でことなるか、検定を行いたいのですが、なかなか方法が見つからず調べていたところ、青木先生のホームページにたどり付きました。
以下記載いたしました質問 (1), (2) ご教授いただけましたら幸いです。

下の表は対象生物170個体をクラスター解析し、
各群(A〜F)で各クラスター(V1〜V9)に分配された個体数です。
	V1	V2	V3	V4	V5	V6	V7	V8	V9
A群 6 6 4 3 4 9 3 8 6
B群 2 1 8 4 0 1 1 2 3
C群 9 2 1 1 1 0 0 4 1
D群 8 2 3 2 0 3 4 2 3
E群 1 0 1 0 1 0 0 0 0
F群 3 3 4 9 6 4 4 7 10
質問(1)「クラスターに含まれる個体数の傾向は群間で異なるか」検定を行いたいのですが、
その場合、カイ二乗検定→ライアンの方法で多重比較 を行う方法は正しいでしょうか?

また、青木先生の説明およびRコードを基に、カイ二乗検定で全体としての差が認められたことを確認後、ライアン法をひとまず行ってみたのですが、
エラーメッセージ「Error: n > 0 は TRUE ではありません 」と表示され、検定を行うことができませんでした。
質問(2)ライアン法は0が含まれていると行えないのでしょうか?

(1), (2) についてご解答のほどよろしくお願いいたします。

Re: クラスタリングされた個体の割合の比較とライアンの方法について
  投稿者:青木繁伸 2019/06/10(Mon) 09:47 No. 22740
ライアン法は,あなたのデータには使えません(使用例参照)

全ての2群の組み合わせ(A:B, A:C, ..., E:F の15組〕について 2 x 9 の行列で chisq.test を行い,p 値をBonferroni法で解釈する(0.05/15)

手作業でやるのは大変なので,以下のようにすればよい
x = matrix(c(
6, 6, 4, 3, 4, 9, 3, 8, 6,
2, 1, 8, 4, 0, 1, 1, 2, 3,
9, 2, 1, 1, 1, 0, 0, 4, 1,
8, 2, 3, 2, 0, 3, 4, 2, 3,
1, 0, 1, 0, 1, 0, 0, 0, 0,
3, 3, 4, 9, 6, 4, 4, 7, 10),
nrow = 6, byrow=TRUE)

g = function(ij) {
i = ij[1]
j = ij[2]
y= rbind(x[i,], x[j,])
z = colSums(y)
y2 = y[, z != 0]
p = chisq.test(y2)$p.value
cat(i, j, p, "\n")
return(c(i, j, p))
}

gn = nrow(x)
cn = ncol(x)
p0 = 0.05 / choose(gn, 2)
res = data.frame(t(combn(6, 2, g)))
colnames(res) = c("i", "j", "p.value")
res$sig = ifelse(res$p.value > p0, "n.s.", "*")
print(res)
結果は以下のようになる
   i j     p.value  sig
1 1 2 0.060309971 n.s.
2 1 3 0.095874639 n.s.
3 1 4 0.365146072 n.s.
4 1 5 0.558144611 n.s.
5 1 6 0.384082770 n.s.
6 2 3 0.034933166 n.s.
7 2 4 0.234596686 n.s.
8 2 5 0.246627056 n.s.
9 2 6 0.177712096 n.s.
10 3 4 0.267205617 n.s.
11 3 5 0.438290867 n.s.
12 3 6 0.005751997 n.s.
13 4 5 0.160816574 n.s.
14 4 6 0.077873765 n.s.
15 5 6 0.441102843 n.s.

Re: クラスタリングされた個体の割合の比較とライアンの方法について
  投稿者:PV 2019/06/10(Mon) 12:05 No. 22741
青木先生、

早速のご解答、ご丁寧な解説を頂き誠にありがとうございました。
Rで回して無事回答を得られることができました。

引き続き質問となり大変恐縮ですが、ライアン法について
もう少しお伺いさせていただいてもよろしいでしょうか?

今回私のデータではライアン法は適していないとのことでしたが、
私のデータは、郡ごとに各クラスターにあてはめられた個体数の頻度をみており、
青木先生のライアン法およびチューキー法の説明にある
「いくつかの群の比率に全体として差が見られたときに,どの群間に差があるかを検定する手法。(比率の差の多重比較(対比較))」
にあてはまると(素人的に)思ってしまったのですが
私の解析は「比率の差」を見ているわけではないということでしょうか?

また、私のデータは以下のホームページにある例題(ここから青木先生のHPを知りました)

猫も杓子も構造化
http://nekomosyakushimo.hatenablog.com/entry/2017/08/26/222728

のカテゴリー数を増やしただけだと思ったのですが、この解釈は間違っているということでしょうか?
統計解析が素人のため愚問で大変恐縮です。

今後解析を行う上でぜひ参考にさせていただきたいためご解答いただけましたら幸いです。
どうぞよろしくお願いいたします。

Re: クラスタリングされた個体の割合の比較とライアンの方法について
  投稿者:青木繁伸 2019/06/10(Mon) 12:15 No. 22742
ライアン法,チューキー法は,群(グループ)の数を k としたとき,k x 2 の行列データを対象にするのです。例えば,正否とか失敗・成功のような2値データが対象です。一方(もう一方)の「比率」が対象になります。

一般的に k x m 行列は,群により「分布が違うか」を対象にすることになります。

「カテゴリー数を増やしただけ」ではありません。

Re: クラスタリングされた個体の割合の比較とライアンの方法について
  投稿者:PV 2019/06/10(Mon) 18:42 No. 22743
青木先生、

なるほど、大変勉強になりました。
お忙しいところ大変丁寧にご説明頂き誠にありがとうございました。


R textでのイタリック表示
  投稿者:コロン 2019/05/26(Sun) 17:39 No. 22733
お世話になっております。

R plotで描画をしておりますが,グラフ内に p<0.01と表示させたいと思っております。

text(1.5, 70, labels="p <- 0.01)

pだけをイタリックで表示させたいのですが,どのようにすればよろしいでしょうか。

ご教示頂けますでしょうか。

Re: R textでのイタリック表示
  投稿者:青木繁伸 2019/05/26(Sun) 21:08 No. 22734
実行可能な最小限のテストプログラムを添付してください

Re: R textでのイタリック表示
  投稿者:コロン 2019/05/27(Mon) 08:38 No. 22735
青木先生

申し訳ございませんでした。簡単なコードではございますが,よろしくお願いいたします。

m1 <- 30
m2 <- 40
x <- c(1,2)
allm <- c(m1, m2)
plot(x, allm, ylim=c(0, 100))
text(1.5, 60, labels="p < 0.05")

Re: R textでのイタリック表示
  投稿者:青木繁伸 2019/05/27(Mon) 18:13 No. 22736
できばえはちょっといまいちなのですが(字が薄い),以下のようにすれば可能
m1 <- 30
m2 <- 40
x <- c(1,2)
allm <- c(m1, m2)
plot(x, allm, ylim=c(0, 100))
text(1.5, 60, labels="p", vfont = c("sans serif", "italic"))
x2 = 1.58 # この値を locator で決めるか,試行錯誤で決める
text(x2, 60, labels=" < 0.05", vfont = c("sans serif", "plain"))
または,書き始めの位置指定を文字列の左端として(pos=4),2 番目の文字列の最初を p のための空白を一つ確保し,
text(1.5, 60, labels="p", vfont = c("sans serif", "italic"), pos=4)
text(1.5, 60, labels=" < 0.05", vfont = c("sans serif", "plain"), pos=4)
とする(この方が試行錯誤不要なのでよいかも)。
"sans serif" の代わりに "serif" も試して,良いと思う方を使ってください。

Re: R textでのイタリック表示
  投稿者:コロン 2019/05/27(Mon) 20:52 No. 22737
青木先生

お忙しい中,ありがとうございました。今作成しておりますグラフにも使用することができました。

[1] [2] [3] [4] [5] [6] [7] [8]

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