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

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

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

２項分布の発見者
　　投稿者：統計の子猫 2020/05/09(Sat) 17:30 No. 22911
２項分布を初めて見つけた人と時代を探しております。色々と調べているのですが、見つけることが出来ず、ここでお聞きすることにしました。

ポアソン分布であれば、1838年にシメオン・ドニ・ポアソンが提案したということが分かっていますし、正規分布であれば、実はガウスよりも前にフランスの数学者ド・モアブルが見つけたと言われたりしています。
ド・モアブルは1730年代に、２項分布のnを大きくしていくと、２項分布の形が正規分布で近似できることを発見したと書いてあります。このことから、２項分布の発見は正規分布よりも古いだろうと予想しますが、誰がどのようにして２項分布を発見したり気付いたのかが書かれた文献などを見つけることが出来ません。
どなたか、２項分布の発見者または提案者が誰なのか、いつごろ発見されたのか、お教えいただけないでしょうか。また、２項分布の歴史が分かる文献なども教えていただけると助かります。
よろしくお願いいたします。
Re: ２項分布の発見者
　　投稿者：青木繁伸 2020/05/09(Sat) 22:13 No. 22912 HomePage
ヤコブ・ベルヌーイ　？
Wikipedia の記事
Re: ２項分布の発見者
　　投稿者：統計の子猫 2020/05/10(Sun) 09:36 No. 22913
青木先生、

コメント、有難うございます。ご紹介いただいた記事を読みました。
ベルヌーイが、ベルヌーイ試行をまとめた場合、つまり二項分布を考えなかったとは思えません。二項分布は、1713年頃にベルヌーイによって発見されたと言えそうです。
有難うございました。


【Ｒ】相関係数行列から距離行列への変換式
　　投稿者：明石 2020/05/01(Fri) 13:44 No. 22906
青木先生　様；

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

掲示板において、たしか、
相関係数行列（類似度）から距離行列（非類似度）への変換式について、掲載されていた記憶があり、
ワード検索しながら、過去分も調べていますが、見つかりません。

Google先生にお聞きすると、
距離=2*(1-相関係数)
というサイトが見つかりますが、
掲示板で、青木先生からご説明のあった変換式とは異なる、という印象をもっています。

青木先生が掲示板でご紹介されていた、変換式を使いたいと思います。

相関係数行列（類似度）から距離行列（非類似度）への変換式について、
この記事No.だよ、と教えていただければ、大変に助かります。

毎々、お手数をおかけいたします。
//
Re: 【Ｒ】相関係数行列から距離行列への変換式
　　投稿者：青木繁伸 2020/05/02(Sat) 21:49 No. 22909 HomePage
http://aoki2.si.gunma-u.ac.jp/lecture/misc/clustan.html
でしょうか？

相関係数行列から平方距離行列を求める場合について書いています（距離行列ではありません）。

> 変数のクラスター分析を行う場合には，変数 i と変数 j の相関係数を rij

> としたとき，2 変数間の距離が次式で表されることになるので，個体のクラスター分析と同じように取り扱うことができる。

d2ij=2 (1−rij),　　(i,j=1,2,…,p)
御礼（Re: 【Ｒ】相関係数行列から距離行列への変換式）
　　投稿者：明石 2020/05/03(Sun) 08:47 No. 22910
青木先生　様；

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

今回も、助けていただきました。

ありがとうございます。
距離行列という、私の記憶違いでしたが、
私が、探していたのは、この記事です。

ありがとうございました。
お手数をおかけしてしまいました。
//


Rの出力について
　　投稿者：もも 2020/05/01(Fri) 15:21 No. 22907
Rの出力方法について，基本的なことで恐縮ですが，ご教授いただけたら幸甚です。

a <- 1
b <- 3
c <- 0.5
min(a,b,c)

とすると，「0.5」と出力されますが，これを「c」と出力するにはどのようにすれば良いでしょうか。
数日間あれこれと調べてみましたが，どうしても分からず，質問させていただきました。
よろしくお願いします。
Re: Rの出力について
　　投稿者：青木繁伸 2020/05/02(Sat) 21:25 No. 22908 HomePage
何をしたいのか本質がよく理解されないのですが。
min_value_name = function(a, b, c) {
	par = c(deparse1(substitute(a)),
		deparse1(substitute(b)),
		deparse1(substitute(c)))
	x = min(a,b,c)
	print(par[which(x == c(a, b, c))])
}
のような関数を定義すると，
> a <- 1
> b <- 3
> c <- 0.5
> min_value_name(a, b, c)
[1] "c"
> x <- 2
> y <- 1
> z <- 5
> min_value_name(x, y, z)
[1] "y"
のような結果を得ることは出来ます，しかし
> min_value_name(3, 2, 1)
[1] "1"
のような使い方をすると，望む結果は得られません。

もし，引数の名前が必ず a, b, c であるなら，以下のような関数を定義することも出来ます。
f <- function(a, b, c) {
	names = letters[1:3]
	return(names[which(c(a, b, c) == min(c(a, b, c)))])
}
これだと，以下のような使用例になります
> f(1,2,3)
[1] "a"
> f(20, 10, 40)
[1] "b"
> f(100, 200, 10)
[1] "c"



【Ｒ】 apply系関数を使った、文字列の連結
　　投稿者：明石 2020/04/29(Wed) 11:07 No. 22903
青木先生　様；

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

標題について、ご教示をいただけましたら、大変に助かります。
毎々、お手数をおかけいたします。
何卒、どうぞよろしくお願いいたします。

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

2本のベクトル x y があります。

x <- c("りんご", "みかん", "桃", "バナナ", "西瓜")
y <- c(0, 2, 0, 2, 1)

要素の個数は同じです。

ベクトル要素の番号どうしに、対応関係があります。

x の各文字列を、y の数値だけ繰り返して、連結した文字列を作成したいと思います。

イメージ的には、以下のようになります。

"りんご"*0, "みかん"*2, "桃"*0, "バナナ"*2, "西瓜"*1

結果として、以下の文字列ベクトルを得たいと思います。

"みかん" "みかん" "バナナ" "バナナ" "西瓜"

以下、コードを作成しました。

s <- NULL
for (i in seq_along(x)) {
v <- rep(x[i], y[i])
s <- c(s, v)
}

これを、apply系関数で書き直すとすれば、どのようになるのか、
ご教示をいただけましたら、大変に助かります。

毎々、お手数をおかけいたします。
何卒どうぞ、よろしくお願いします。
//
Re: 【Ｒ】 apply系関数を使った、文字列の連結
　　投稿者：青木繁伸 2020/04/29(Wed) 21:14 No. 22904 HomePage
複数の引数を取れる mapply を使うと良いと思います。
unlist と unname を使う必要があります。
unname(unlist(mapply(function(a, b) rep(a, b), x, y)))
御礼（Re: 【Ｒ】 apply系関数を使った、文字列の連結）
　　投稿者：明石 2020/04/30(Thu) 10:21 No. 22905
青木先生　様；

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

今回も、新しい知識を教えてくださり、誠にありがとうございました。

mapply()関数の名前は聞いたことがありますが、　
自分が悩んでいた問題で、mapply()関数の使い方をご教示いただけましたので、
よく理解できました。

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


【Ｒ】リストからデータフレームへの変換（その２）
　　投稿者：明石 2020/04/27(Mon) 17:26 No. 22900
青木先生　様；

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

先にご教示いただきました、 No. 22896 【Ｒ】リストからデータフレームへの変換
について、
追加のご教示をいただけましたら、大変に助かります。

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

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

以下の例題について、
apply系関数を使う、リストからデータフレームへの変換方法を、
ご教示をいただけましたら、大変に助かります。

s <- c("りんご_apple_林檎", "みかん_orange_蜜柑", "バナナ_banana", "メロン_melon")
s <- strsplit(s, "_")

リストに格納された4本のベクトルの要素の個数は、異なります。

> s[[1]]
[1] "りんご" "apple" "林檎"
> s[[2]]
[1] "みかん" "orange" "蜜柑"
> s[[3]]
[1] "バナナ" "banana"
> s[[4]]
[1] "メロン" "melon"

やりたいことを、以下にお示しをします。

> s[[1]][1]
[1] "りんご"
> s[[2]][1]
[1] "みかん"
> s[[3]][1]
[1] "バナナ"
> s[[4]][1]
[1] "メロン"
のように、s[[*]][1]のベクトル要素だけを集めて、ベクトルに格納する。

同様に、
s[[*]][2]のベクトル要素だけを集めて、ベクトルに格納する。
s[[*]][3]のベクトル要素だけを集めて、ベクトルに格納する。

以下のような3本のベクトルを得ることができましたら、データフレームを作成できます。
> x1
[1] "りんご" "みかん" "バナナ" "メロン"
> x2
[1] "apple" "orange" "banana" "melon"
> x3
[1] "林檎" "蜜柑" NA NA

s[[*]][1]のベクトル要素だけを集めて、ベクトルに格納する。
s[[*]][2]のベクトル要素だけを集めて、ベクトルに格納する。
s[[*]][3]のベクトル要素だけを集めて、ベクトルに格納する。
の操作で、もし、apply系関数が利用できるならば、
その方法をご教示いただけましたら助かります。

毎々、お手数をおかけいたします。
何卒、どうぞよろしくお願いいたします。
//
Re: 【Ｒ】リストからデータフレームへの変換（その２）
　　投稿者：青木繁伸 2020/04/27(Mon) 22:08 No. 22901 HomePage
以下のようにすれば，もっともらしいとは思いますが，スマートでもなんでもないと思います。
素直にプログラミングする方がよいかと思います。
s <- c("りんご_apple_林檎", "みかん_orange_蜜柑", "バナナ_banana", "メロン_melon")
v  <- sapply(s, function(u) {
	str <- character(3)
	a <- unlist(strsplit(u, "_"))
	len <- length(a)
	str[1:len] = a
	str})
df <- data.frame(t(v))
rownames(df) <- NULL
df[df==""] <- NA
df
df[,1]
df[,2]
df[,3]
実行結果
> df
      X1     X2   X3
1 りんご  apple 林檎
2 みかん orange 蜜柑
3 バナナ banana 
4 メロン  melon 
> df[,1]
[1] "りんご" "みかん" "バナナ" "メロン"
> df[,2]
[1] "apple"  "orange" "banana" "melon" 
> df[,3]
[1] "林檎" "蜜柑" NA     NA
御礼（Re: 【Ｒ】リストからデータフレームへの変換（その２））
　　投稿者：明石 2020/04/28(Tue) 09:12 No. 22902
青木先生　様；

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

今回も、ご丁寧で、ありがたいご教示をいただき、とてもうれしいです。
心から、心より、御礼を申し上げます。
ありがとうございました。

以下、補足です。
------------------------------------

s[[*]][1]のベクトル要素だけを集めて、ベクトルに格納する、
という個所を、　例えば、sapply(s, FUN,1) 　などのように、

s[[*]][2]のベクトル要素だけを集めて、ベクトルに格納する、
という個所を、　例えば、sapply(s, FUN,2) 　などのように、

もし、上記のように書けたとすると（書けるかどうか、その知識がありません）　
FUNは、どのように書けるのか、という強い興味がございました。

投稿の動機が不純だったような気がしており、反省しています。
やはり、プログラミングも、素直さが大切です。
//


【Ｒ】リストからデータフレームへの変換
　　投稿者：明石 2020/04/25(Sat) 17:09 No. 22896
青木先生　様；

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

リストからデータフレームへの変換について、ご教示をいただけましたら、大変に助かります。

-------------------------------
以下、簡単な例でお示しします。

s <- c("りんご_apple", "みかん_orange", "バナナ_banana", "メロン_melon")
s <- strsplit(s, "_")

> s
[[1]]
[1] "りんご" "apple"

[[2]]
[1] "みかん" "orange"

[[3]]
[1] "バナナ" "banana"

[[4]]
[1] "メロン" "melon"

上記リストから、以下の2列からなるデータフレームを作成したいと思います。
lhs rhs
1 りんご apple
2 みかん orange
3 バナナ banana
4 メロン melon

上記リストを、unlist()関数でベクトルに変換して、
ベクトル要素の奇数番目と偶数番目を、
作成したいデータフレームの変数 lhs と rhs に代入することでも得られますが、

ご教示いただきたいことは、リストからデータフレームへの直接変換の方法です。

お手数をおかけいたします。
どうぞよろしくお願いいたします。
//
Re: 【Ｒ】リストからデータフレームへの変換
　　投稿者：青木繁伸 2020/04/25(Sat) 17:44 No. 22897 HomePage
「直接変換」という意味がよく分かりません。
df <- data.frame(t(matrix(unlist(s), 2)))
ではだめなんですか？
御礼（Re: 【Ｒ】リストからデータフレームへの変換）
　　投稿者：明石 2020/04/25(Sat) 19:08 No. 22899
青木先生　様；

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

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

「直接変換」という、紛らわしい表現で、ご迷惑をおかけいたしました。
大変に失礼をいたしました。

青木先生からご教示いただきました方法で、問題解決できました。
今回も、大変に良い勉強をさせていただきました。

私が「直接変換」という表現を使った期待感として、
apply, map, reduce などの関数を使う方法があるとすれば、どのように書くのだろう、
という思いがありました。

「直接変換」という、紛らわしい表現で、ご迷惑をおかけいたしました。
大変に失礼をいたしました。
//


【Ｒ】文字列ベクトルの突き合わせ
　　投稿者：明石 2020/04/14(Tue) 12:14 No. 22892
青木先生　様；

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

文字列ベクトルの突き合わせについて、ご教示をいただけましたら、大変に助かります。

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

以下の3本の文字列ベクトルがあります。

v <- c("リンゴ", "ミカン", "バナナ", "スイカ", "メロン", "モモ")
s1 <- c("バナナ","モモ")
s2 <- c("リンゴ","メロン")

s1、s2を、vと突き合わせをして、以下の行列を作成したいと思います。

【出来上がり】
リンゴ ミカン バナナ スイカ メロン モモ
s1 0 0 1 0 0 1
s2 1 0 0 0 1 0

> charmatch(v, s1)
[1] NA NA 1 NA NA 2

> charmatch(v, s2)
[1] 1 NA NA NA 2 NA

この結果を、ifelse()関数で、0/1に置換したいと思い、色々とやっていますが、うまくいきません。

ご教示をいただけましたら助かります。
お手数をおかけいたします。
//
Re: 【Ｒ】文字列ベクトルの突き合わせ
　　投稿者：青木繁伸 2020/04/14(Tue) 14:21 No. 22893 HomePage 
v <- c("リンゴ", "ミカン", "バナナ", "スイカ", "メロン", "モモ")
s1 <- c("バナナ","モモ")
s2 <- c("リンゴ","メロン")
#
a1 <- charmatch(v, s1) # NA NA 1 NA NA 2
ifelse(is.na(a1), 0, 1) # 0 0 1 0 0 1
# または
1 - is.na(a1) # 0 0 1 0 0 1
#
a2 <- charmatch(v, s2) # 1 NA NA NA 2 NA
ifelse(is.na(a2), 0, 1) # 1 0 0 0 1 0
# または
1 - is.na(a2) # 1 0 0 0 1 0
# 以上を踏まえ
x <- 1 - rbind(is.na(a1), is.na(a2))
rownames(x) <- c("s1", "s2")
colnames(x) <- v
x
#    リンゴ ミカン バナナ スイカ メロン モモ
# s1      0      0      1      0      0    1
# s2      1      0      0      0      1    0
または %in% を使って
y = rbind(v %in% s1, v %in% s2) + 0
rownames(y) <- c("s1", "s2")
colnames(y) <- v
y
御礼（Re: 【Ｒ】文字列ベクトルの突き合わせ）
　　投稿者：明石 2020/04/14(Tue) 15:11 No. 22894
青木先生　様；

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

今回も、大変に良い勉強をさせていただきました。
心から、心より、御礼を申し上げます。
ありがとうございました。
//


【Ｒ】文字列ベクトルの突き合わせの高速化
　　投稿者：明石 2020/04/11(Sat) 13:44 No. 22889
青木先生　様；

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

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

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

2つの文字列ベクトルの突き合わせの高速化について、
ご教示をいただけましたら、大変に助かります。

使うデータは、2種類です。
・文字列ベクトル v
・参照テーブル tbl

・文字列ベクトル v は、
　自然言語処理で、テキストから切り出した単語の集合であり、長大です。
　表記ゆれが多くありますので、統制の必要性があります。

　表記ゆれの例示としては、"りんご","リンゴ","林檎","アップル","apple"など

・参照テーブル tbl は、
　上記の表記ゆれを統制するために、参照する表です。

　tbl は2列からなります。
　1列目は、参照文字列
　2列目は、置換文字列

　tbl の例示としては、
　りんご 林檎
　アップル 林檎
　apple 林檎
　など

やりたいことは、以下です。
・文字列ベクトル v を、for()関数でループを回す。
・v の各要素を、参照テーブル tbl の1列目「参照文字列」と照合する。
　ここで、標題の「文字列ベクトルの突き合わせ」が出てきます。
　== 演算子　を使いました。
・もし、v の要素が、tbl の1列目「参照文字列」にあれば、
　そのv の要素を、 tbl の2列目「置換文字列」で置換する。

for()関数の2重ループが原因と思いますが、
気が遠くなるほどの長い処理時間となりますので、
文字列ベクトルの突き合わせの高速化について、
ご教示をいただけましたら、大変に助かります。

お手数をおかけいたします。
//
Re: 【Ｒ】文字列ベクトルの突き合わせの高速化
　　投稿者：青木繁伸 2020/04/11(Sat) 21:34 No. 22890 HomePage
以下のようなプログラムだと，実行時間はいかほどになりますか？

dict = matrix(c(
"りんご", "林檎",
"アップル", "林檎",
"apple", "林檎",
"レモン", "檸檬",
"れもん", "檸檬",
"lemon", "檸檬",
"remon", "檸檬"), byrow=TRUE, nc=2)
v = c("りんご", "アップル", "apple",
	"檸檬", "レモン", "れもん", "lemon", "remon",
	"イチゴ", "スイカ")
a = unname(sapply(v, function(s) {
	i = which(dict[, 1] == s)
	ifelse(length(i) == 0, s, dict[i, 2])
}))
実行結果
> a
 [1] "林檎"   "林檎"   "林檎"   "檸檬"   "檸檬"  
 [6] "檸檬"   "檸檬"   "檸檬"   "イチゴ" "スイカ"
御礼（Re: 【Ｒ】文字列ベクトルの突き合わせの高速化）
　　投稿者：明石 2020/04/12(Sun) 09:53 No. 22891
青木先生　様；

お忙しいところを失礼いたします、明石と申します。
今回も、大変にありがたい、貴重なご教示をいただきました。
ありがとうございました。

劇的な改善をはかることができました。

エラップス時間の比較です。
・青木先生の方法：16秒
・私の従来方法　：25分経過しても、まだ終わりません..........

データ量が多くなれば、この効果は極めてありがたいです。

青木先生からご教示いただきました貴重なノウハウを、勉強しています。

心から、心より、御礼を申し上げます。
//


【Ｒ】文字列のコード判定（ひらがな、漢字）
　　投稿者：明石 2020/04/03(Fri) 16:11 No. 22886
青木先生　様；

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

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

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

3つの文字列があります。

x <- c("りんご","林檎","apple")

1番目の要素は、ひらがな　のみ
2番目の要素は、漢字　　　のみ

どのようにしたら、
1番目の要素は、ひらがな　のみ
2番目の要素は、漢字　　　のみ
という判定ができるのでしょうか。

お手数をおかけいたします。
何卒どうぞよろしくお願いいたします。
Re: 【Ｒ】文字列のコード判定（ひらがな、漢字）
　　投稿者：青木繁伸 2020/04/03(Fri) 18:21 No. 22887 HomePage
Python での記事ですが，
https://note.nkmk.me/python-re-regex-character-type/
が参考になると思います。

R では以下のようになります。
HIRAGANA = "[\u3041-\u309F]" # ひらがなの正規表現
grepl(HIRAGANA, "あ") # TRUE
grepl(HIRAGANA, "ア") # FALSE

KATAKANA = "[\u30A0-\u30FF]" # カタカナの正規表現
grepl(KATAKANA, "ア") # TRUE
grepl(KATAKANA, "あ") # FALSE

# 簡易
KANJI0 = "[\u4E00-\u9FFF]" # 漢字の正規表現（簡易版）
grepl(KANJI0, "漢") # TRUE
grepl(KANJI0, "々") # FALSE
grepl(KANJI0, "〆") # FALSE
grepl(KANJI0, "〇") # FALSE

# 完全
KANJI = "[\u2E80-\u2FDF\u3005-\u3007\u3400-\u4DBF\u4E00-\u9FFF\uF900-\uFAFF\U00020000-\U0002EBEF]"  # 漢字の正規表現（完全版）
grepl(KANJI, "漢") # TRUE
grepl(KANJI, "々") # TRUE
grepl(KANJI, "〆") # TRUE
grepl(KANJI, "〇") # TRUE

# 複数文字の場合，すべてがひらがな，カタカナ，漢字かどうかの判定
is.all.hiragana = function(s) {
  t = grepl(HIRAGANA, strsplit(s, "")[[1]])
  # print(t)
  all(t)
}

is.all.katakana = function(s) {
  t = grepl(KATAKANA, strsplit(s, "")[[1]])
  # print(t)
  all(t)
}

is.all.kanji = function(s) {
  t = grepl(KANJI, strsplit(s, "")[[1]])
  # print(t)
  all(t)
}

is.all.hiragana("ひらがな")
is.all.hiragana("ひらがなイガイ")

is.all.katakana("カタカナ")
is.all.katakana("カタカナいがい")

is.all.kanji("漢字")
is.all.kanji("漢字いがい")
御礼（Re: 【Ｒ】文字列のコード判定（ひらがな、漢字））
　　投稿者：明石 2020/04/03(Fri) 19:58 No. 22888
青木先生　様；

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

今回も、ご丁寧で、分かりやすい、ご教示をいただき、誠にありがとうございました。

私はPythonも勉強し始めましたので、
青木先生がご紹介してくださいましたページも、あわせて、参照いたします。

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


共変量で調整した二元配置分散分析
　　投稿者：田崎 2020/02/05(Wed) 08:57 No. 22885
お世話になっております。
恐れ入りますが、以下の共変量で調整した二元配置分散分析について、検定手法が妥当であるかご教示いただけますと幸いです。

独立変数：BMI（3群）と飲酒状況（3群）
従属変数：血液の生化学値
共変量：年齢

二元配置分散分析を行った結果、BMIと飲酒それぞれの独立変数はP<0.001、両者の相互作用はP＝0.75となりました。
そこで、BMI、飲酒状況それぞれで、3群間の差を解析したいと考えました。
この場合、SPSSでは共変量を設定したままではpost-hoc testができないため、BMIと飲酒それぞれについて共分散分析（多重比較としてBonferroni法を採用）をpost-hoc testとして行うのが良いのではと考えました。

お時間がありますときにご教示いただけますと幸いです。
どうぞよろしくお願いいたします。


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))
}
両方の言語を知っていれば，移植もそう難しくはないのでは？


Ｒ　区間分割処理
　　投稿者：明石 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の２つを引数として渡すだけで、
findInterval(x, c(-0.03000, 10.68895, 37.48632, 64.28368))
を自動的に組み立てて、実行して、所属する階級を得ることができたら、大変に助かります。

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

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

ご教示をいただけましたら大変に助かります。
新年早々、お手数をおかけいたします。
//
Re: Ｒ　区間分割処理
　　投稿者：aoki 2020/01/05(Sun) 12:34 No. 22880
findInterval(x, mx[,1]) だけでよいのでは？
Re: Ｒ　区間分割処理
　　投稿者：明石 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
> 自分の結果の解釈が正しいのか

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


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

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

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

----------------------------------------
時系列データで、Ｒに組み込みの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: Ｒ　時系列データの変換
　　投稿者：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: Ｒ　時系列データの変換）
　　投稿者：明石 2019/12/28(Sat) 18:26 No. 22875
青木先生　様；

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

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

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

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

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


サーストンの一対比較法で求めた間隔尺度から有意差の有無を調べる方法
　　投稿者：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
初めて投稿させていただきます。お手数をおかけしますが、ご教授いただけましたら幸いです。
質問紙調査において、

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

という設問（内容は一部改変しています）で、設問２と４で「はい」「いいえ」の割合が異なるかを調べたいと思っています。
この場合、設問２と４の回答者の中には、「設問２だけ答えた人」「設問４だけ答えた人」「設問２と４ともに答えた人」が混在していることになるかと思います。
この場合、設問２と４のデータは対応があるデータなのか、対応がないデータなのか（設問２と４ともに答えた人は対応があると思いますが）、どのように考えたらよいのでしょうか。
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
お世話になっております。ご教授ください。

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

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

よろしくお願いいたします。
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
教育心理分野で、ある実習後の学習意欲を目的変数とし、実習中の環境要因２つ、実習生の個人要因２つを説明変数として重回帰分析をしようと考えていました。

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

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. 目的変数、説明変数とも正規分布している必要がある。

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

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

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


計算関係にあるデータ間の相関検定の可否
　　投稿者：コダマ 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   


※答えを求めているわけでは断じてありません。ご理解下さい。
以下のデータが得られたとして　行う分析とそれにより得られる結論を述べよというものです。

会社Ａ 資本金1000万 従業員100名
業績平均値：総務部７０ 営業部５０ 開発部５０

会社Ｂ 資本金3億 従業員400名
業績平均値：総務６０ 営業３０ 開発６０
いずれの平均値に対してもＳＤ＝5とする

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

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







Re: どう扱えば良いか困っています


　　投稿者：青木繁伸  2019/08/25(Sun) 17:25  No. 22815   


業績平均値は各部の部員の業績の平均値ということでしょうかね？
資本金の違いは関係ないのでは？
全従業員数の違いも関係ないと思いますが，総務部，営業部，開発部 の人数がわからないので，分析不可能では？
つまり，会社全体の業績平均値を知ろうにも，各部の人数がわからなければ計算しようがない。

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







Re: どう扱えば良いか困っています


　　投稿者：りゅうへい  2019/08/25(Sun) 18:36  No. 22818   


ご回答ありがとうございます。
問題には　仮想の業績得点の平均値と書かれていました。
平均値の認識は　ご指摘の通りで合っていると思います。
3部署合計人数＝従業員数なのか　3部署以外の職員を含め　従業員数なのかが記載されていませんでした。　
3部署合計が全従業員数　という前提であれば　各部の人数が分からずとも
会社全体の平均値は　Ａは　(70+50+50)/3=56.67 B (60+30+60)/3=50　
と　求めることは可能でしょうか？

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

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







Re: どう扱えば良いか困っています


　　投稿者：青木繁伸  2019/08/25(Sun) 20:40  No. 22819   


> 3部署合計が全従業員数　という前提であれば　各部の人数が分からずとも会社全体の平均値は　Ａは　(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が遅くなりお詫び申し上げます。
先生からいただいた説明を何度も何度も読み直してみましたが
私の説明が悪くて質問の意図が伝えられていないと思いました

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

そして同様に上記と独立して、
「4つある形態の一つであるCHという形態を呈する患者が、Ａ型では全体の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 
> 「たとえばＡ型という疾患においてCR, CH型の頻度を比較をする」ということではなく...CRという形態はA型がB型よりも有意に多いことを示せるか

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

> 「4つある形態の一つであるCRという形態を呈する患者が、Ａ型では全体の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)をやる」、というのは
下のような２ｘ２の表を作成して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　の３つのシートからなる
"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   


青木先生　様；

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

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

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











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





- J o y f u l 　N o t e -

Modified by　i s s o