![[RjpWiki]](image/../Rlogo.png "[RjpWiki]")
#author("2024-04-15T06:25:23+09:00","","") #author("2024-04-15T06:37:05+09:00","","") COLOR(green){SIZE(20){初心者のための R および RjpWiki に関する質問コーナー}} COLOR(red){SIZE(16){以下の 4 つの「●項目」のどれかをクリック！！}}~ ---- COLOR(blue){● [[新規投稿欄:>#sinki]]} --- その前に，「[[投稿における注意事項]]」を読んでから！&br; COLOR(blue){● [[このページの目次>#mokuji]]} --- 質問への回答・コメントの参照&br; COLOR(blue){● [[最新のスレッド>#saisin]]} --- 最も最近に書き込まれた「親」記事（最新のコメントが付いた記事ではありません）&br; COLOR(blue){● [[Ｑ＆Ａ (初級者コース) 過去の目次]]} --- 今までの 17 個の書庫の参照&br; ---- &aname(mokuji);COLOR(green){●● 目次 ●●} 参照は，個々の項目をクリック（質問への回答・コメントの参照も，個々の項目をクリック）~ #contents ---- &aname(sinki); COLOR(green){●● 新規投稿 ●●} まずは，「[[投稿における注意事項]]」を読んで，それに従ってください。~ COLOR(red){これに従わない場合には，回答が得られにくい場合があります。} #article **カイ二乗検定について [#pdeea0ee] >[[カイ]] (2024-04-14 (日) 14:17:43)~ ~ よろしくお願いいたします。~ カイ二乗検定を行い、P値を取り出すとき以下のよう、マトリックスにして計算しております。~ ~ d <- c(10, 5, 8, 7)~ d <- matrix(d,nrow=2,byrow=T)~ d <-chisq.test(d)~ d$p.value~ ~ しかし、~ ~ data <- matrix(c(10, 5, 8, 7, 16, 23, 21, 8, 27, 44, 61, 22, 10, 15, 20, 24, 38, 76, 23, 77), nrow=5, byrow=T)~ colnames(data)<-c("a","b","c","d")~ ~ のようなデータがある場合、一気にカイ二乗検定を行い、ｄ列の横にｐ値を入れていきたいのですが、マトリックスが出来なく困っております。~ どの様にすればよいのかご教示いただけると幸いです。~ （R version 4.2.2, win11)~ // - 質問者様のやろうとしていることに対して、ご提示のような変数dataに元データを入れておくのが適しているか、という点は一旦おいておくことにしますが……各行の4つ組の値を2x2の表とみなして、それに対するp値を最後の列に追加したい、ということなら下記のようにできますね。前述のとおり、このデータとp値の組を無理やりこのような5列の変数にすること自体の是非は正直ありますが。素直にやるなら、データは2x2xNのarrayにしておいて、1:Nまでの検定の結果のp値は別変数に長さNのベクトルとしておく、というのがstraight forwardでしょうか。 -- &new{2024-04-15 (月) 06:25:23}; - 質問者様のやろうとしていることに対して、ご提示のような変数dataに元データを入れておくのが適しているか、という点は一旦おいておくことにして……各行の4つ組の値を2x2の表とみなして、それに対するp値を最後の列に追加したい、ということでよいでしょうか？それでしたら下記のように可能です。~ > print(data) a b c d [1,] 10 5 8 7 [2,] 16 23 21 8 [3,] 27 44 61 22 [4,] 10 15 20 24 [5,] 38 76 23 77 > p <- apply(data, 1, function(x) { chisq.test(matrix(x, nrow = 2))$p.value }) > cbind(data, pval = p) a b c d pval [1,] 10 5 8 7 7.093881e-01 [2,] 16 23 21 8 2.012252e-02 [3,] 27 44 61 22 1.955183e-05 [4,] 10 15 20 24 8.518854e-01 [5,] 38 76 23 77 1.287840e-01 前述のとおり、2x2のデータ群と、そのそれぞれに対する検定結果のp値を、無理やりこのような5列の行列にすること自体の是非は正直あります。（一覧性の都合など、質問者様の固有の理由があるのだと思いますが。）~ プログラム的に素直にやるなら、データは2x2xNのarrayにしておいて、1:Nまでの検定の結果のp値は別変数に長さNのベクトルとして得る、というのがstraight forwardでしょうか。 -- &new{2024-04-15 (月) 06:25:23}; > data2 <- array(t(data), dim = c(2, 2, nrow(data))) > print(data2) , , 1 [,1] [,2] [1,] 10 8 [2,] 5 7 , , 2 [,1] [,2] [1,] 16 21 [2,] 23 8 , , 3 [,1] [,2] [1,] 27 61 [2,] 44 22 , , 4 [,1] [,2] [1,] 10 20 [2,] 15 24 , , 5 [,1] [,2] [1,] 38 23 [2,] 76 77 > apply(data2, 3, function(x) { chisq.test(x)$p.value }) [1] 7.093881e-01 2.012252e-02 1.955183e-05 8.518854e-01 1.287840e-01 #comment **library(imager) library(magick)の間で画像をやりとりしたい。 [#p183510f] >[[たたた]] (2023-11-20 (月) 13:23:54)~ ~ お世話になっております。表記のとおり２つのパッケージ間で画像をやりとりしたいです。画像をimagerで用いるcimgで読み込み、magickへの読みこむことは以下の方法でできるのですが、逆にmagickの画像からcimgにする方法をどなたか教えていただければ幸いです。~ image1 <- imager::load.image("image1.png") image1.mag <- magick::image_read(as.cimg(image1)) 一応、手元で、以下のようにすると、それっぽくはなるのですが、値が数値ではなくrawだし、次元数も変わってしまっており、どうしたらよいものかと思っております。~ image_matrix <- magick::image_data(image1.mag) image_imager <- imager::cimg(image_matrix) どうぞよろしくお願いします。~ // #comment **with関数でformulaオブジェクトを代入するとエラーが出る [#d05e0b06] >[[徒弟]] (2023-08-28 (月) 21:26:18)~ ~ 皆様のお知恵をお貸しください~ 以下のコマンドの通りを試したところ、(3)でのみエラーが出ます~ ~ -------------------~ ~ set.seed(1234)~ tmp <- data.frame(Outcome=c(rep(1,10),rep(0,10)),X=rnorm(20))~ tar_formula <- as.formula(paste0("Outcome","~X"))~ ~ #(1)~ with(tmp,glm(Outcome~X,family=binomial))~ ~ #(2)~ with(tmp,glm(as.formula(paste0("Outcome","~X")),family=binomial))~ ~ #(3)~ with(tmp,glm(tar_formula,family=binomial))~ ~ -------------------~ ~ エラーメッセージは" eval(predvars, data, env) でエラー: ~ オブジェクト 'Outcome' がありません"です ~ どういよう理由によってエラーが発生するのか教示頂ければ幸いです~ // - さらにコメント申し上げます。本当の狙いはmiceなどによって生成したデータを上記tmpとして、(3)を実施したいと思っていましたが、同様のエラーが発生しました。調査したところmiceとか関係なく、上記の状況で同様のエラーが出るため、この点を解決できないかと考え相談させていただきました。 -- [[徒弟]] &new{2023-08-28 (月) 21:57:48}; - さらにコメント申し上げます。本当の狙いはmiceなどによって生成したデータを上記tmpとして、(3)を実施したいと思っていましたが、同様のエラーが発生しました。調査したところmiceとか関係なく、上記の状況で同様のエラーが出るため、この点を解決できないかと考え相談させていただきました。 -- [[徒弟]] &new{2023-08-28 (月) 22:22:39}; - as.formulaのenv引数がparent.frame()ですので、事前に実行しておいたtar_formula <- as.formula(hoge)は、R_GlobalEnvを参照しにいっていて、withの張った空間内のOutcomeを見つけられない、ということではないでしょうか？ > tar_formula <- as.formula("Outcome ~ X") > str(tar_formula) Class 'formula' language Outcome ~ X ..- attr(*, ".Environment")=<environment: R_GlobalEnv> with内でas.formulaを実行した場合には、つくられたformulaはwithがつくった空間内を参照しにいっています。 -- &new{2023-08-29 (火) 06:39:21}; > with(data.frame(), str(as.formula("Outcome ~ X"))) Class 'formula' language Outcome ~ X ..- attr(*, ".Environment")=<environment: 0xa6bb540> - 早速ご教示いただき誠にありがとうございます。よくわかりました。 -- [[徒弟]] &new{2023-08-29 (火) 12:50:06}; - Rで?formulaを実行し、Environmentsのセクションを見たことでより詳細に状況を確認できました。しかしなんで参照先を決める必要があるのでしょうね？ -- [[徒弟]] &new{2023-08-29 (火) 13:22:55}; - 質問者様の想定は「Outcomeっていう変数名ということだけ指定しておいて、その実体は実際の実行時に柔軟に（現環境から）読み取ってよ」ということだと思います。が、裏を返すとそのような「柔軟」な挙動は、そのformulaオブジェクトを実際にevalする時点でのEnvironmentに依存して、予想外の別の変数を読みにいってしまう危険がある、ともいえるかと思います。なので、そういう間違えが起こらないよう、「formulaをつくる時点で、その計算対象の変数が普通はもうあるのだから、それがどのEnvironmentにあるかも指定しておく」という設計思想なのでしょうか。~ 別の言い方をすると、何らかの「計算の結果」を返すのが一般的な関数だとしたら、as.formulaという関数は「これから計算すべきこと」を返す関数なので、そういう配慮が必要といえるかもしれません。もっと端的に、プログラムの仕様上Environmentが一意に指定されていないとマズい挙動をする例を挙げられそうな気もしますが、Rの内部実装について詳しくないので他の詳しいかたにお任せします。 -- &new{2023-08-29 (火) 13:45:41}; - 私の想定していることを具体化いただきありがとうございます。仰る通りです -- [[徒弟]] &new{2023-08-29 (火) 21:08:31}; - 私の想定していることを具体化いただきありがとうございます。仰る通りです -- [[徒弟]] &new{2023-08-29 (火) 22:33:25}; #comment **ベイジアンネットワーク解析dealの矢印について [#x3ed22ae] >[[晴猫]] (2023-08-27 (日) 19:36:31)~ ~ yokkunさんのページ~ https://yokkuns.hatenadiary.org/entry/20110928/1317164851~ にならってベイジアンネットワーク分析を試みたところ、ネットワークの矢印が正しく表示されません。~ mybanlist <- matrix(c(~ 2,1, 3,1, 4,1, 5,1, 6,1, 7,1, 3,2, 4,2, 5,2, 6,2, 7,2, 3,4, 4,3), ncol = 2, byrow = TRUE) banlist(ht.nw) <- mybanlist~ で矢印を禁止するノードを指定します。その後~ plot(ht.nw)~ でネットワークを描いてみると、矢印がないと指定したノードの間に赤い破線で矢印が描かれます。しかし、その後の計算は正しくできるようです。ですからこれは矢印の表示の問題だと思います。通常は黒い実線でノード間の因果関係が表示されます。赤い破線は逆方向の因果関係を表しているようです。正方向の因果関係を表示させるにはどうしたらよいのでしょうか。~ OSはWindows11, Rはversion 4.2.1を使用しています。よろしくお願いします~ // - R 4.3.1 on Ubuntu環境でやってみましたが、再現しませんねぇ……引用元のサイト様と同じ解析結果となります。質問者様の環境での実行の標準出力と描画結果などを貼り付けてみてはいかがでしょうか？~ また、この解析についてはまったく知らないのであくまで説明文を読んでの想像ですが、この禁止ノードの指定というものには方向性があるのですよね？質問者様の環境でどのような描画になっているのか不明ですが、その「赤い破線」というのが逆方向の因果を表わすのであれば、banlistで指定した向きでは禁止できてない、といったことはないのでしょうか？（たとえばbanlistで両方向の関係性をBANすれば消える、とか。）~ 依然として、そもそも引用元サイト様と違う解析結果になっていること自体の原因は不明ですが…… -- &new{2023-08-28 (月) 06:01:55}; - 実行の出力です -- [[晴猫]] &new{2023-08-30 (水) 00:44:02}; - ご返信ありがとうございます。状況が再現して、質問者様の意図がわかりました。（おかしなページに結果が貼られてしまっているのは、ご自身でどうにかしていただくこととして……）~ 前提として、わたしは依然このネットワーク分析というものについて無知ですので、その点はご了承ください。さて、このライブラリのautosearch関数は、実行中に分析結果を段階的に描画していき、最後に残った画面が最終分析結果となるようです。ここで質問者様は「plot(ht.nw)をし直すと赤い矢印が引かれてしまう」と仰っています。（これ自体は再現します。）しかしこのht.nwという変数はautosearchする前のinitial networkなので、autosearchした結果は持っていません。問題の赤い破線の矢印は、banlistとして指定したエッジそのものです。「ht.nwには、これらのエッジは想定しないものとして指定してありますよ」ということを意味しています。そのエッジが分析の結果残った、ということではありません。~ autosearchの結果は、今回のサンプルコードだとht.searchという別のオブジェクトに保存していますよね？このオブジェクトのnwという名前の要素に、分析後のネットワークの情報が保存されているようです。 > names(ht.search) [1] "nw" "table" "trylist" > class(ht.search$nw) [1] "network" ですので、plot(ht.nw)ではなくplot(ht.search$nw)とすれば、autosearchし直さなくても解析結果を再表示できます。plot(ht.nw)は分析前のデータの描画、plot(ht.search$nw)が分析後の結果の描画です。その際、件の赤矢印が「事前のbanlistの内容」として一緒に描画されますが、それが邪魔ならplot時のshowban引数をFALSEにすれば描画しないことが可能です。 plot(ht.search$nw, showban = FALSE) networkクラスのオブジェクトの描画には、plot関数に用意されたnetworkメソッドが使われています。plot.networkのヘルプを参照していただければ、描画の仕様はわかるかと思います。 -- &new{2023-08-30 (水) 07:26:23}; - ご回答ありがとうございます -- [[晴猫]] &new{2023-08-30 (水) 21:12:34}; - ご回答ありがとうございます -- [[晴猫]] &new{2023-08-30 (水) 21:14:32}; - ご回答ありがとうございます -- [[晴猫]] &new{2023-08-30 (水) 21:16:54}; #comment **conda-forge の r-base パッケージの win-64 のバージョンについて [#kd2e8ddd] >[[136s]] (2023-08-02 (水) 11:10:23)~ ~ 質問: 掲題のバージョンが現時点で v4.1.3 ですが、linux-64 や osx-64 は v4.3.1 になっています。win-64 がアップデートされないのは、何か理由があるのでしょうか？~ 技術的不可能性があるのか、需要がないのか、、など。~ ~ 参考: https://anaconda.org/conda-forge/r-base~ ~ 経緯: Windows 10 マシンで、conda 環境内で R を管理して Python から rpy2 で R を使おうとしています。~ 通常、`conda install -c conda-forge r` コマンドで conda 環境内に R をインストールすると思います。~ ただ、現在の最新の [[r パッケージ:https://anaconda.org/conda-forge/r]] が noarch で v4.3 であり、r-base のバージョン不一致になってしまいます。~ これ自体は r パッケージのインストール時に v4.1 を指定すれば解決するのですが、そもそも r-base が v4.1.3 からバージョンアップされないのか気になったので質問を投稿いたしました。~ // - そもそも、conda-forge の r, r-base を誰が管理してるのかもよくわからず…。Google Search で conda site:www.r-project.org と検索してもそれらしい情報はなさそうです。 -- [[136s]] &new{2023-08-02 (水) 11:15:33}; - 自己解決しました。 https://github.com/conda-forge/r-base-feedstock/issues/248 から関連情報を確認していくと、CRAN が 2022-04 から C ライブラリとして URCT を使用することで、conda-forge との整合が付かなくなる？可能性があるようで、2022-03 リリースの v4.1.3 で止まってるみたいでした。 -- [[136s]] &new{2023-08-02 (水) 12:05:52}; #comment **RStudio のパッケージits.analyasisの解析方法について [#k8a63535] >[[masa]] (2023-06-04 (日) 14:03:04)~ ~ ①結果の入力方法(他に必要なパッケージなど) ②結果の解釈 を知りたいです。最低限①だけでも助かります。 // #comment **ggplot2による積み上げ棒グラフの各棒に合計値を入れたい [#p24f5645] >[[島田泰夫]] (2023-06-02 (金) 17:16:52)~ ~ 教えてください~ ~ ggplot2で積み上げ棒グラフを作図したとき、棒の上に合計値を表示したいのです~ ~ tableのクロス集計表（re）で~ ~ 0 1 2 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 402 0 3 2 1 2 37 0 0 0 0 0 1 10 15 20 0 0 1 0 2 1 0 0 404 0 0 0 0 0 0 1 2 2 12 1 5 3 2 0 1 9 0 0 0 0 0 0 406 0 2 0 0 15 1 19 21 3 0 3 8 0 0 0 0 0 2 0 0 0 0 0 411 0 0 0 20 5 6 10 11 12 8 9 0 0 0 2 6 10 3 2 0 2 1 1 416 16 11 18 1 0 0 1 0 2 0 0 0 0 0 0 8 0 0 14 16 1 1 0 418 0 0 0 0 0 0 0 20 13 5 2 2 10 18 62 94 60 11 3 4 4 27 0 422 0 0 0 0 0 0 0 0 9 6 5 10 12 14 8 30 26 22 14 9 3 4 3 428 0 0 15 4 9 1 12 2 3 2 2 15 1 10 2 4 0 0 0 0 9 16 6 429 10 6 0 0 0 21 7 2 7 0 0 1 17 16 9 4 3 0 3 1 4 0 3 が得られたので、データフレーム変換して ~ re2 <- as.data.frame(re) Var1 Var2 Freq 1 402 0 0 2 404 0 0 3 406 0 0 .... ~ となりましたので~ ~ g <- ggplot(re2, aes(x = Var2, y = Freq, fill = Var1)) g <- g + geom_bar(stat = "identity") g <- g + geom_text(aes(label = Freq, y = Freq + 3), size = 2) plot(g) ~ とすると、積み上げる「各要素」の数値が入ってしまいます~ これを避けて、積み上げた「合計値」を各棒の上部に記載したいのですが~ どのようにすればよろしいでしょうか？~ // 自己解決できそうです re3 <- as_tibble(re2) re4 <- re3 %>% group_by(Var2) %>% summarise( total= sum(Freq)) ggplot(data = re4, aes(x = Var2, y = total)) + geom_bar(stat = "identity", fill = "grey", colour = "black") + geom_text(aes(label = total, y = total + 3), size = 2) + geom_bar(data = re2, aes(x = Var2, y = Freq, fill = Var1), stat = "identity") -- [[島田泰夫]] &new{2023-06-02 (金) 18:30:05}; #comment **PCAの結果の分析について [#o31c462c] >[[R初心者]] (2023-05-08 (月) 09:58:23)~ ~ 特徴量が1000個あるデータに対してPCAを実行しました．有効な特徴量を抽出するには，どの計算結果を参考にすればよろしいでしょうか．．．基本的な質問で申し訳ありません．~ // - Rでのプログラミングについてお答えすれば解決できる、というご質問ではないように感じます。そもそもPCAがなにを算出していてどのように使われる解析方法なのか、といった部分は統計の教科書等で充分にご理解されていてのご質問でしょうか？もしそうであれば、具体的にコードを実行するうえで何が疑問となっているのか、もうちょっと詳しくお書きになられたほうが回答が得られやすいかと思います。 -- &new{2023-05-08 (月) 11:33:25}; - コメントいただきありがとうございます． 「Rでまなぶデータサイエンス２多次元データ解析法」で勉強しながらRのコーディングをするとともに で出力結果が何を意味しているのか調べています 75×945行列のデータに対して2クラスに分けるために有効な特徴量を945個の中から探しています． 累積寄与率が99％を超えるのは第５３主成分（全主成分は７５）まで必要になっているとともに， 主成分得点と主成分負荷量をから有効な特徴量を選択するためのデータの解釈法が分かりません． -- [[R初心者]] &new{2023-05-08 (月) 12:32:03}; - 主成分をもとに「有効な（元データの）変量を選ぶ」というのは、ちょっとよくわからないのですが...。元データのもつ共変動関係をできるかぎり低次元に縮約するというのが主成分分析の根幹だと思うので、そのうえで「どれかひとつ（ないしいくつか）の特徴量に立ち戻る」という感覚がわかりません。また、第53主成分まで取るというのがそもそもわたしの知っている主成分分析の感覚とだいぶ違う感じです。（結局、人間の解釈可能な次元数に縮約しきれていないということなので。）分野によるのかもしれません。質問者様の分野でそうした解析がメジャーであるなら、そのような解析を使っている先行研究の例示などは可能でしょうか？2クラスにわけてそれぞれの特徴を知りたいということなら、クラスタリングでもして、各変量のグループごとの分布などを比較するというのがストレートかなと、個人的には思ってしまいますが。 -- &new{2023-05-08 (月) 12:47:05}; - コメントいただきありがとうございます．私も第３主成分で99％近くの累積寄与率になるかと思っていたのですが，60％程度の累積寄与率でした．もう少し主成分分析について勉強いたします．お手数ですが主成分分析や次元削減についておすすめの参考書などありましたら，ご教示いただけますと幸いです． -- [[R初心者]] &new{2023-05-08 (月) 12:54:25}; - 本日コメントをいただきまして，私の中で，一番不明な点は，次元削減を行った主成分を，次の学習にどのように利用するのかが分かっていないことが分かりました． -- [[R初心者]] &new{2023-05-08 (月) 13:46:02}; - 本日コメントをいただきまして，私の中で，一番不明な点は，次元削減を行った主成分を，次の学習にどのように利用するのかが分かっていないことが分かりました． -- [[R初心者]] &new{2023-05-08 (月) 13:51:43}; - うーんと……「次の学習」というのが、機械学習みたいなことを意味してらっしゃるのか、あるいは単に「次の段階としてのデータの解釈」的なことをおっしゃっているのか、またわたしには分かりかねますが…。まぁそれはおいておくとして。~ 一般論として元データの変量間に共変動関係がなければ次元は縮約しようがないわけですから、累積寄与率があまりあがっていかないということは、次元縮約向きではないデータなのかもしれません。変量が945もあるのに、レコードが75しかないというのも、単純に厳しそうにも思えます。（サンプルをクラス分けしたいというだけなら、やっぱりクラスタリングのほうが向いているのではないでしょうか？変量の内容にもよりますが。）~ とはいえ、前述のとおりこれらは分野によるでしょう。わたし自身はばいおろじーの人間ですが、Social Science系だとまた事情が違うだろうと思います。やはり、なさろうとしている解析と似た先行研究を探して、そこで使われている手法を参考にするのが早そうです。主成分分析自体の勉強については、データマイニング系の成書にはたいがい載っているし、何でされてもよいかと思います。 -- &new{2023-05-11 (木) 07:03:35}; - コメントいただきありがとうございます．「次の段階としてのデータの解釈」という意味での学習です．言葉足らずで申し訳ありません．主成分分析の基礎から再勉強しようと考えております．また，クラスタリングについても検討いたします。今後ともよろしくお願いいたします． -- [[R初心者]] &new{2023-05-11 (木) 15:57:05}; - 機械学習屋さんなのでしょうか。「2クラスに分ける」というのなら，（変数選択をする）判別分析，ロジスティック判別分析とかをやるべきではないでしょうか?[75×945行列]はデータ規模としては異常。サンプルサイズ75は小さすぎるし，変数個数945は異常に大きすぎる。-- &new{2023-05-16 (火) 23:14:39}; - わたしも（自分にとっては）馴染みのない感じのデータだなぁとは思いつつお返事していました。が、分野によって「ありとあらゆる変数を取りまくること自体は自動的に簡単にできる」みたいな事情があるとか、サンプルサイズにしても本データ集め前のパイロットデータであるとか、まあ特有の事情などがあったりするのかなということで。テスト的なデータを使って、Rを勉強がてら、いろんな解析を試したりされてるところなのかな？とか。 -- &new{2023-05-17 (水) 07:21:29}; #comment **Rcmdrのインストールとエラーについて [#wee008f5] >[[hs0519]] (2023-03-14 (火) 23:18:04)~ ~ 具体的な質問というより、一般的な考え方の質問です。RのWindows版をインストールした直後、Rを起動させて、install.packages("Rcmdr",dependencies=TRUE)などとやってインストールするのですが（ダウンロード先はTokyo）、エラーが出ます。インストール作業のプロセスが進行していることは動作をみていてわかりますね。メッセージなどが出て進行状況がわかります。そして最後にエラーが出るという図式です。~ この流れの中で異質なこと、独特なところは何もありません。エラーが出ないときもあります。使用者側はここまでの流れで工夫するところがないので、提供側の問題としか思えません（ダウンロード中の欠損もありうる？）。そのためエラーメッセージを見ながらフィックスする努力は徒労でかつ、生産的でないように感じられます。やはりエラーメッセージをたよりにフィックスするのか、Ｒそのものの更新を待つべきか悩みが出ます。この部分は果報は寝て待てなんじゃないかと思うのですが、実際どうなのでしょうか。~ // - （わたしも具体的回答ではなく個人的感想ですが）普段からLinuxを使っていると、そういう「エラー出たから足りないライブラリ足さなきゃ」的な作業に慣れ（毒され？）すぎていて、あまり疑問を抱いたことがありませんでした。よく考えたら、そりゃもちろん徒労感ですよね。「エラーが出ないときもある」というのが「再現性がなくエラーの有無が変わる」という意味だとすれば、なおさら。単に「depedenciesを全自動で解消できるパッケージとできないパッケージがある」という意味であれば、（良い悪いはともかく）まぁありうるかなと思いますが。 -- &new{2023-03-16 (木) 16:09:28}; - かつて(20年程前)の Code set Independence ( Win9xやEUC ) 等に基づいたコードが今も生きていて, そのせいだと思います. ユーザー名がASCIIコードのみであれば動作するのではないでしょうか? -- [[なかま]] &new{2023-03-23 (木) 21:08:01}; - Rまわりの最初期から活躍されてる方のコメントの心強さときたら。。。 :) そんなところに原因が潜んでいる可能性があるものなのですね。 -- &new{2023-03-24 (金) 07:07:00}; #comment **mFDというパッケージのエラーについて [#k8d38d29] >[[まつ]] (2023-03-14 (火) 20:56:57)~ ~ mFDというパッケージでmFD::quality.fspacesを行った際に~ ~ Error in if (maxdim_pcoa < 1) { : the condition has length > 1~ ~ というエラーが出ました。教えて頂けると嬉しいです。ちなみにmaxdim_pcoaは２以上の色々な値で試していますが同じエラーが出ます。~ // - 直前の質問への回答と被りますが、「投稿における注意事項」にしたがい、minimal reproducible exampleをご提示いただけませんか？パッケージ自体がおかしそうか、使いかたがおかしそうかは、とりあえずエラーを再現してみることができるとあたりを付けやすいのです。 -- &new{2023-03-16 (木) 16:11:50}; #comment **指数の作り方 [#j19a4a36] > (2023-03-13 (月) 14:40:09)~ ~ 日付、性別、年代、検索行動１、検索行動２、検索行動３、収益額(この値は99%が0で、そのほかは1000円から50000円くらいの値が入っています。)~ という情報があったときに収益をもたらすかどうかを示す指数を作りたいです。~ ~ 指数の条件としては以下の2つです。~ ・値が大きいほど収益額1000円超の顧客の割合は大きくなる傾向にある~ ・値が大きいほど収益額が大きくなる傾向にある~ ~ lightgbmのregressionでやってみても全然だめでした。~ お手数ですがお知恵をお貸しください。~ // - 「投稿における注意事項」にあるとおり、何をどうやったらどうダメだったのかが再現できるようなコードやデータ（minimal reproducible example）をご提示いただかないと、お力になれないかと思います。現状では、その「指数」というもののつくりかたに問題があって「ダメ」だったのか、コードにエラーがあって「ダメ」だったのかすらわかりません。もっといえば、元のデータにそもそも収益額の予測能力がなければ、収益額を予測する「指数」が実際つくれない（「ダメ」という結果になるのが正しい）のかもしれません。 -- &new{2023-03-14 (火) 07:02:47}; #comment **Rからffmpegを呼び出す方法 [#i3f1da84] >[[hs0519]] (2023-03-10 (金) 08:02:35)~ ~ saveVideo({~ for (i in 1:5){ plot(runif(10), ylim = c(0,1)) } }, interval = 1.0, video.name = "TEST.html")~ というスクリプトだけで、ffmpegが呼び出せるとのことなので、実行してみると、~ The command ""ffmpeg"" is not available in your system. Please install FFmpeg~ というエラーメッセージが出ます。~ ffmpegはインストールされており、パスも通っているのでシェルのコマンドラインから起動できます。~ Rからffmpegを起動するにはどのようにしたらいいでしょうか。~ ~ バージョンにセンシティブな問題でしょうか。~ Rは4.2.2, ffmpegは5.1.2が入っています。ただし、ffmpegはいろんなところ(irfanviewとか)から入る可能性があるので5.1.2でない可能性が大です。~ パスが先になっている方だと思いますが。~ // - saveVideo()のヘルプ内の例では、ani.options("ffmpeg")で指定されるのは.exeまで含んだFFmpegへのフルパスになっています。システムのエイリアスが認識されていないのが問題であれば、とりあえずフルパスで実行ファイルを指定してみてはどうでしょう？ちなみにわたしは、シェルのコマンドをpaste()でつくってからsystem()で実行してffmpegを使っています。OS依存などをある程度回避できるので、saveVideo()のような既存関数を使うより汎用性が高いように思います。 -- &new{2023-03-10 (金) 12:35:08}; - ありがとうございます。フルパスは保存先フォルダ名だけでなく、実行ファイル名まで含んで指定するわけですね。 -- [[hs0519]] &new{2023-03-10 (金) 14:31:04}; - ヘルプではそのような記法になっているので、とりあえずそれで試してみるのがストレートかな、と。なお「保存先」のディレクトリはani.options("imgdir")に相対パス指定、と指示されていますね。 -- &new{2023-03-10 (金) 15:51:30}; #comment **UTFのRでWindowsの文字をよみたい。 [#c7c6cbb4] >[[たたた]] (2023-03-01 (水) 06:05:33)~ ~ エクセルからコピーしたタブ区切りの「データ０ データ１（ほげ） データ２？ データ３? データ４」というデータがあり、 name <- read.table("clipboard",sep="\t")という処理を行いたいと思います。 手元のRの環境で、（や？があっても読めるものと読めない物があるのですが、読めないほうでも読むにはどのようにしたらよいかわかりますでしょうか。 読めるものは：R version 3.6.2 (2019-12-12)、LocaleがJapanese_Japan.932 。読めないほうは：R version 4.2.1 (2022-06-23 ucrt)。LocaleがJapanese_Japan.utf8 。いずれも、 Windows 10 x64 (build 19041)です。 なお、文字コードをShift-JIS系に指定しても以下のエラーが出てしまいます。 read.table("clipboard",sep="\t",fileEncoding = "CP932") 警告メッセージ: 1: file(file, "rt", encoding = fileEncoding) で: 引数 'encoding' は無視されます 2: read.table("clipboard", sep = "\t", fileEncoding = "CP932") で: incomplete final line found by readTableHeader on 'clipboard' ~ // - read.table("clipboard", sep = "\t", fileEncoding = "CP932") で: incomplete final line found by readTableHeader on 'clipboard' これは，エンコーディング等とは無関係なエラー。入力行の最後に改行コードがありませんということ。多くの場合はちゃんと読めているけど，クリップボードにコピーするときに改行を含んでコピーすればよい。 -- &new{2023-03-04 (土) 10:23:52}; #comment **as.*() と as_*() 関数の違いについて [#p4924198] >[[OKD]] (2023-02-05 (日) 15:27:12)~ ~ こんにちは．~ 初歩的な質問で大変恐縮ですが，以下のご回答をよろしくお願い致します．~ ~ 早速ですが， as.*() と as_*() の違いが理解できません．~ "*" には，matrixやvector，factor などを想定しています．~ "as"の部分は，”is”の場合でも違いがわかっていません．~ どなたかご教示いただけますでしょうか？~ // - Rのビルトインの関数はas.matrixやas.factorであり、個人的にはas_matrixといった名前は見たことがないです。ご利用されている特定のパッケージ内で定義されている関数、といった可能性はありませんか？ -- &new{2023-02-06 (月) 06:37:44}; - 早速，ご回答ありがとうございます．ご指摘の通りかもしれません．corrrやpurrrなどのパッケージ特有の関数ということなんですね．ネット上で紹介されているコードを写すにしても，見かけ上，as.matrixやas.vectorと同じ挙動？をしていたので，何故ビルトイン関数を用いないのかとても不思議でした． -- [[OKD]] &new{2023-02-06 (月) 21:08:37}; - corrrパッケージのas_matrixの挙動を調べてみました。同関数のヘルプに書いてあるとおり「correlate関数で作成された『cor_dfクラス』のデータをmatrixに変換する」というはたらきのようですね。 > x <- correlate(iris) Non-numeric variables removed from input: `Species` Correlation computed with • Method: 'pearson' • Missing treated using: 'pairwise.complete.obs' > class(x) # cor_dfという独自のクラスとなる [1] "cor_df" "tbl_df" "tbl" "data.frame" > print(x) # printするとmatrixぽく表示されるが… # A tibble: 4 × 5 term Sepal.Length Sepal.Width Petal.Length Petal.Width <chr> <dbl> <dbl> <dbl> <dbl> 1 Sepal.Length NA -0.118 0.872 0.818 2 Sepal.Width -0.118 NA -0.428 -0.366 3 Petal.Length 0.872 -0.428 NA 0.963 4 Petal.Width 0.818 -0.366 0.963 NA > str(x) # 実際はtermという文字列要素をもつdata.frame cor_df [4 × 5] (S3: cor_df/tbl_df/tbl/data.frame) $ term : chr [1:4] "Sepal.Length" "Sepal.Width" "Petal.Length" "Petal.Width" $ Sepal.Length: num [1:4] NA -0.118 0.872 0.818 $ Sepal.Width : num [1:4] -0.118 NA -0.428 -0.366 $ Petal.Length: num [1:4] 0.872 -0.428 NA 0.963 $ Petal.Width : num [1:4] 0.818 -0.366 0.963 NA > as_matrix(x) # 数値matrixに変換 Sepal.Length Sepal.Width Petal.Length Petal.Width Sepal.Length NA -0.1175698 0.8717538 0.8179411 Sepal.Width -0.1175698 NA -0.4284401 -0.3661259 Petal.Length 0.8717538 -0.4284401 NA 0.9628654 Petal.Width 0.8179411 -0.3661259 0.9628654 NA このcor_dfクラスというオブジェクトは、基本的にはdata.frameなのですが、変量名にあたる情報をtermという要素としてもっています。なので、単にビルトインのas.matrixを適用すると、文字列型のmatrixに強制変換されてしまうようです。 > as.matrix(x) # 4x5の文字列matrixになってしまった term Sepal.Length Sepal.Width Petal.Length Petal.Width [1,] "Sepal.Length" NA "-0.1175698" " 0.8717538" " 0.8179411" [2,] "Sepal.Width" "-0.1175698" NA "-0.4284401" "-0.3661259" [3,] "Petal.Length" " 0.8717538" "-0.4284401" NA " 0.9628654" [4,] "Petal.Width" " 0.8179411" "-0.3661259" " 0.9628654" NA それを防ぐため、「cor_dfオブジェクトをmatrixへと、数値データを数値のまま変換する」ための独自関数がas_matrixである、という感じでしょうか？~ こういう場合は、as.matrixにcor_dfクラス用のメソッドを追加するというのが定石かと思っていたのですが……最近のRのコーディングのお作法がどのようになっているかあんまりフォローできてないので、よくわかりません。でも「cor_dfという独自クラス専用で、as.matrixのはたらきをする関数」に、as_matrixという汎用っぽい名前がついているのは、やっぱりちょっとミスリーディングな気もします。（じっさい質問者さんが混乱されたように。）~ ちなみにas_matrixは「cor_dfクラス以外には普通のas.matrixのはたらきをする」わけではなく、あくまでcor_dfクラス用になっているみたいですね。 > y <- data.frame(a = 1:3, b = 4:6) > print(y) a b 1 1 4 2 2 5 3 3 6 > as.matrix(y) a b [1,] 1 4 [2,] 2 5 [3,] 3 6 > as_matrix(y) UseMethod("as_matrix") でエラー: 'as_matrix' をクラス "data.frame" のオブジェクトに適用できるようなメソッドがありません 混同して使ってしまうとまずいです。 -- &new{2023-02-07 (火) 07:38:32}; - 検証までしてくださり，本当にありがとうございます．違いまで，よく分かりました．危うくどっちでも同じなら．．．的なノリで使ってしまいそうでした． -- [[OKD]] &new{2023-02-07 (火) 19:34:22}; #comment **当日0:00:00からの経過時間 [#s1d6d3c5] >[[asao]] (2022-12-16 (金) 14:26:56)~ ~ こんにちは~ 素人質問で恐縮ですが時系列データの取り扱いについての質問です。~ ~ 現在以下の形式のようなデータを用いて、時刻と使用率の関係を推定しようとしています。~ ~ ー－－－－－データー－－－－－－ー－－－－－－－－－－ー－－~ beginning_time ocp_rate inflow outflow~ 2020/6/1 23:58:00 0.023255858 1 0~ 2020/6/1 23:59:00 0.893255895 0 1~ 2020/6/2 0:00:00 0.023255826 1 2~ 2020/6/2 0:01:00 0.103255844 0 0~ 2020/6/2 0:02:00 0.045255836 0 1~ 2020/6/2 0:03:00 0.113255814 2 0~ ー－－－－－データー－－－－－－－－－－－－－－－－－－－－~ そこで新たに~ その日付の深夜0時（1行目でいうと2020/6/1 00:00:00）~ からの経過時間（分）を示す行を追加したいと考えております。~ ~ どのように計算すればよいでしょうか。~ ~ 今のところ、各日付の00:00:00を示す行を追加して、beginning_time列から引き算しようと試みておりますが、各日付の00:00:00の列を作ることにも苦戦しております。~ ご教授いただけますと幸いです。~ // - 経過時間（分）の「列」、ということでよいでしょうか？いろいろなやりかたがあると思いますが……時間をあつかう場合には、POSIXlt型またはPOSIXct型を使うのがとりあえずは定石かなと思います。~ もとのデータテーブルの仕様がわからないので、とりあえず文字列matrixをもとに作業してみると、こんなカンジでしょうか。 dat <- matrix(c( "2020/6/1 23:58:00", "0.023255858", "1", "0", "2020/6/1 23:59:00", "0.893255895", "0", "1", "2020/6/2 0:00:00", "0.023255826", "1", "2", "2020/6/2 0:01:00", "0.103255844", "0", "0", "2020/6/2 0:02:00", "0.045255836", "0", "1", "2020/6/2 0:03:00", "0.113255814", "2", "0"), byrow = TRUE, ncol = 4) colnames(dat) <- c("beginning_time", "ocp_rate", "inflow", "outflow") time <- as.POSIXlt(dat[, "beginning_time"], format = "%Y/%m/%d %H:%M:%S") h <- as.numeric(format(time, "%H")) m <- as.numeric(format(time, "%M")) elapsed <- h * 60 + m dat <- cbind(dat, elapsed = elapsed) print(dat) 今回ご質問の作業「だけ」をするなら、sub関数などで時:分だけ取り出したりする方法もあります。が、時間データをもとに「何月何日の」とか「何曜日のデータだけ」とかアレコレいじる作業をしようと思うと、最初にPOSIXltなどにしてしまっておくのがラクかなと思います。 -- &new{2022-12-16 (金) 16:29:10}; - ご回答、ありがとうございます。 頂いた回答をもとに、解決することができました。 as.numeric(format(time, "%H")) という書き方で、POSIXltの中の時間の部分を実数で取り出すことができるんですね。 勉強になります。 -- [[asao]] &new{2022-12-16 (金) 16:58:50}; - 解決してよかったです。フォーマットの記述様式については、strptime関数のヘルプなどに載っていますので、必要に応じてご参照ください。（あらかじめ知っておくと、いつか使いたい場面が来たときに思い出せるかもしれません。） -- &new{2022-12-16 (金) 17:42:57}; #comment **３つの係数により変化する関数を、３次元空間で、色の濃淡で、大小を表したい [#t2ab284a] >[[okadasan]] (2022-12-13 (火) 15:32:33)~ ~ お世話になっております。~ まだまだ、Ｒは初心者です。~ http://okadasan.com/S_P.R~ と言う関数を作りました。~ ~ p[r0,gamma,a]と言う関数が、~ r0とgammaとaの変化で変わって行っているのがわかると思います。~ ここで、r0をx軸にして、gammaをy軸にして、aをz軸にして、~ ３次元空間で、それら３つの係数の変化で決まる関数であるp[r0,gamma,a]の各点の大きさを、色の濃さで表したいのです。~ そのような事は、Ｒで出来るのでしょうか？~ plot3dで出来るのかなあ、とか、~ 別なグラフ関数を持って来るのかなあ、とか、~ あるいは、その所作をコードの何処に入れるのかも~ わかりません。~ どなたか、お教えいただけるとありがたいのですが。~ // - うーん、言わんとしていることは分かる気がするのですが、Rで云々以前に、原理的に難しくはないでしょうか？3次元空間の全点を塗りつぶすことになると、手前の点の色で奥の点が塗りつぶされてしまいます。それを「ちょうどよい具合に」透過させるのは、「ちょうどよい具合」を決める方法の定義が必要となります。「3次元空間内で、いくつかの基準値で『等高線の3次元版（曲面になります）』のようなものを描く」というほうが、まだ現実的な気がします。ただ、基準値の決め方や、曲面のメッシュの描き方や密度、3次元上での描画方法など、初心者のかたで簡単にできるものではないように（個人的には）思います。 -- &new{2022-12-13 (火) 16:45:54}; - ご提示いただいた関数の意味や内容の可否はわからないので、とりあえず置いておきます（使ってない変数があったり、わたしには理屈がわからない処理があったりするので）。例えばで、おっしゃっている通りの描画をplot3dでしてみると、こんなカンジでしょうか。。。~ library(rgl) r0 <- 1:20 gamma <- 1:36 a <- 1:100 f <- function(x1, x2, x3) { # just for example # not the true "Sound Pressure" function return((x1 - 10)^2 * abs(x2 - 30) + sin(x3 / 10) * 1000) } xs <- expand.grid(r0, gamma, a) p <- f(xs[, 1], xs[, 2], xs[, 3]) n_col <- 100 i <- round((p - min(p)) / diff(range(p)) * (n_col - 1)) + 1 col <- cm.colors(n_col)[i] plot3d(xs[, 1], xs[, 2], xs[, 3], xlab = "r0", ylab = "gamma", zlab = "a", col = col, aspect = FALSE) 関数（f）の定義の部分は、前述のとおりわたしは理解していないので、適当な算術処理に置き換えました。単に「3次元配列の値を空間内にプロットする例」です。やっぱり空間を密に埋めるので、遮蔽の問題をどうにかしないと実用性がなさそうな気がしてしまいます。ざっと傾向を眺めるだけ、なら何かの役に立つものなのでしょうか。。。 -- &new{2022-12-13 (火) 17:31:28}; #comment **offsetを使った混合モデルにlmerTest::difflsmeans()関数を適用するとエラー [#qa8ec833] >[[noname]] (2022-09-04 (日) 17:59:16)~ ~ lmer()関数でoffsetを指定した混合モデルにlmerTest::difflsmeans()関数を適用すると，offsetで指定した変数名のオブジェクトがないとエラーが出ます．~ 以下で再現すると思います．~ ~ library(lmerTest) data("cake", package="lme4") cake2 <- data.frame(cake,os = round(cake$temp/cake$angle,0)) model <- lmer(angle ~ recipe * temp + (1|recipe:replicate) + os, cake2) model2 <- lmer(angle ~ recipe * temp + (1|recipe:replicate) + offset(os), cake2) anova(model) anova(model2) difflsmeans(model) difflsmeans(model2) # オブジェクト 'os'がありませんというエラー ~ ~ sessionInfo() R version 3.6.3 (2020-02-29) Platform: x86_64-apple-darwin15.6.0 (64-bit) Running under: macOS 10.16 Matrix products: default BLAS: /Library/Frameworks/R.framework/Versions/3.6/Resources/lib/libRblas.0.dylib LAPACK: /Library/Frameworks/R.framework/Versions/3.6/Resources/lib/libRlapack.dylib locale: [1] ja_JP.UTF-8/ja_JP.UTF-8/ja_JP.UTF-8/C/ja_JP.UTF-8/ja_JP.UTF-8 attached base packages: [1] stats graphics grDevices utils datasets methods base other attached packages: [1] lmerTest_3.1-3 lme4_1.1-27.1 Matrix_1.2-18 dbscan_1.1-8 colorout_1.2-2 pacman_0.5.1 loaded via a namespace (and not attached): [1] Rcpp_1.0.9 compiler_3.6.3 pillar_1.8.1 nloptr_1.2.2.2 tools_3.6.3 boot_1.3-24 lifecycle_1.0.1 [8] tibble_3.1.8 nlme_3.1-144 gtable_0.3.0 lattice_0.20-38 pkgconfig_2.0.3 rlang_1.0.4 cli_3.3.0 [15] DBI_1.1.1 dplyr_1.0.9 generics_0.1.3 vctrs_0.4.1 grid_3.6.3 tidyselect_1.1.2 glue_1.6.2 [22] R6_2.5.1 fansi_1.0.3 minqa_1.2.4 ggplot2_3.3.6 purrr_0.3.4 magrittr_2.0.3 scales_1.2.1 [29] MASS_7.3-51.5 splines_3.6.3 assertthat_0.2.1 colorspace_2.0-3 numDeriv_2016.8-1.1 utf8_1.2.2 munsell_0.5.0 ~ stackexchangeやマニュアルを調べたものの，そもそもlmer()関数で実際にoffsetを使っている例*すら見つからず，皆様のお知恵を拝借出来れば幸いです． 一番近い例：https://stats.stackexchange.com/questions/88960/lme4-glmer-problems-with-offset [#pd4b1c0b] ~ %%※投稿における注意事項に記載されているユーザ名の求め方が間違っているようです．%% %% http://www.okadajp.org/RWiki/?%E6%8A%95%E7%A8%BF%E3%81%AB%E3%81%8A%E3%81%91%E3%82%8B%E6%B3%A8%E6%84%8F%E4%BA%8B%E9%A0%85 %% 修正しました． // - 開発元にissueを立てました．動きがありましたらご報告いたします． -- [[noname]] &new{2022-09-05 (月) 12:57:41}; - 自己レスです． -- [[noname]] &new{2022-09-06 (火) 20:15:23}; offsetはformuraの中に入れるのではなく，別の引数として指定するものでした． lmerTestパッケージのマニュアルにlmerでのoffset使用例がないことと，lme4::glmerの書き方を流用したこと，formuraのなかに入れてもエラーが出ないこと，何よりも私の確認不足・能力不足，があわさってお騒がせをしてしまいました． 誠に申し訳ございません． library(lmerTest) data("cake", package="lme4") cake2 <- data.frame(cake,os = round(cake$temp/cake$angle,0)) model2 <- lmer(angle ~ recipe * temp + (1|recipe:replicate), offset=os, cake2) anova(model2) difflsmeans(model2) # エラーは出ない - お疲れさまです。こうした問題解決の記録が残ることでいつか誰かの役に立つかもしれませんし、その集積がWikiの役割だと思うので、自己解決のレスに感謝致します。 -- &new{2022-09-08 (木) 10:09:35}; #comment **brm関数でエラーが起きてしまいます [#rcf07b43] >[[花椒報告]] (2022-09-02 (金) 13:10:19)~ ~ result<- brm(~ formula = sales ~ weather + temperature, family = gaussian(), data = d, seed = 1 )~ ~ を実行すると、~ ~ Error in compileCode(f, code, language = language, verbose = verbose) : ~ C:\rtools42\x86_64-w64-mingw32.static.posix\bin/ld.exe: ~ （以下、10行くらいにわたる長いエラー文）~ が出てしまいます。~ ~ library(rstan)~ library(brms)~ は無事できています。~ ~ 環境は~ R version 4.2.1 (2022-06-23 ucrt)~ Platform: x86_64-w64-mingw32/x64 (64-bit)~ Running under: Windows 10 x64 (build 22000)~ ~ です。~ ~ どのように対処すればよいでしょうか？~ // - 追記です。 pkgbuild::has_build_tools(debug = TRUE) の結果はTRUEなので、Rtoolsはインストール済みだと思います -- [[花椒報告]] &new{2022-09-02 (金) 13:11:30}; - 自己解決しました。 一応報告しておきますと、rstanパッケージは、2.21以上だとエラーが起きると報告されているようです。よって、 -- [[花椒報告]] &new{2022-09-02 (金) 14:24:08}; - 2.19.3をインストールしなおすと、うまくいきました -- [[花椒報告]] &new{2022-09-02 (金) 14:24:35}; #comment **rasterパッケージのcropオブジェクトの使い方について [#s9a8597e] >[[てくてく]] (2022-08-23 (火) 16:18:12)~ R version 4.2.1 (2022-06-23 ucrt)~ Platform: x86_64-w64-mingw32/x64 (64-bit)~ Running under: Windows 10 x64 (build 19044)~ ~~ DEMの描画や、切り取りについて勉強中です。試しに、rasterパッケージ内のcropというオブジェクトのヘルプに記載の内容（下記）、~ ~~ r <- raster(nrow=45, ncol=90)~ values(r) <- 1:ncell(r)~ e <- extent(-160, 10, 30, 60)~ rc <- crop(r, e)~ ~ を実行してみたのですが、rcのextentはeで指定した値ではなく、-160, 8, 30, 58になってしまいます。~ crop後のextentの値が、cropの範囲（上記のe）で指定した値にならないのはどうしてでしょうか。ご教授いただけたらと思います。よろしくお願いいたします。~ // #comment **パッケージ自動起動設定の削除方法について [#i063580f] >[[ん]] (2022-08-04 (木) 14:18:29)~ ~ Mac book air (2018, Intel) OS MontereyでRを使用しています。~ EZRが必要になり、頒布元が提示してくださっている手順に従って、R起動時に自動でEZRも起動する設定を行いました。~ ~ options(Rcmdr=list(plugins='RcmdrPlugin.EZR'))~ library(Rcmdr)~ local({~ old <- getOption('defaultPackages')~ options(defaultPackages = c(old, 'Rcmdr'))~ })~ ~ しかしその後Rstudioを使用するようになり、Rstudio起動→Rにコネクト→EZR自動起動のためxquartzが起動→EZR起動→Rstudioがやっと立ち上がる、という状態になってしまいました。~ ~ options関数について調べてみましたが情報が少なく、~ options(defaultPackages = c(getOption('defaultPackages')))~ と入れてみましたが、自動起動設定は変わりませんでした。~ ~ 起動時設定を初期化する方法はありますでしょうか？~ ご教授いただければ幸いです。~ // - 単に.Rprofileファイルを削除するなり、上記の行をすべてコメントアウトするなりではダメなのでしょうか？ -- &new{2022-08-05 (金) 07:45:11}; - もしかしたら.Rprofileの仕組み自体に慣れていらっしゃらない、ということなのかな…？そのEZRという奴のインストールページの指示に従うことで、質問者さんは自身のホームディレクトリに.Rprofileというファイルを作成し、そのなかに上記のコマンド一式を書き込んだことになります。（ただしドットファイルが非表示とかになっていると、Finderでは見えなかったりするかもしれません。）Rは起動時に、この.Rprofileファイルを自動的に読み込んで実行するという仕様をもっているので、その内容が毎回起動時に実行されて、EZRが呼び出されています。ですので、前述のインストールページの最下部に書かれているような方法でこの.Rprofileファイルを呼び出して中身を消したり、あるいはファイルごと削除したりすれば、自動的にEZRが呼び出されることはなくなるかと思います。 -- &new{2022-08-05 (金) 07:57:19}; - お返事いただいたお二方、ありがとうございます。 -- [[ん]] &new{2022-08-05 (金) 09:37:39}; - ご指摘の通りR初心者で、.Rprofileそのものをよく知らずに指示通り打ち込んでしまいました。ターミナルからの操作で削除し、呼び出されなくなったことを確認できました。.Rprofileについて勉強してみます。ありがとうございました。 -- [[ん]] &new{2022-08-05 (金) 09:40:02}; #comment **モンテカルロ検定 [#xb4bcaa3] >[[はえ]] (2022-07-28 (木) 23:27:37)~ ~ 各⾃に割り当てられた、個体の区画ごとの個体数のデータを使い、そのデータが５％⽔準でランダム分布と有意な差があるといえるかを、R のプログラムを作成して使⽤し、モンテカルロ検定する。以下の３つは必ずつける：(1)その検定結果、(2)モンテカルロ検定を⾏う R の⾃作プログラム、(3)モンテカルロ・サンプルの検定統計量のヒストグラム。~ データ~ 2，0,0,5,0,0,0,3,4,18,0,3,2,6,0,4,0,2,0,3,3,9,5,0,0，15，0，6，0，4，0，0，5，0，0，0，0，2，0，9，0，0，0，2，0，0，0，0，0，0，0，5，0，0，5，6，3，0，0，0，0，0，5，0，3，12，5~ 総個体数＝156 区画数＝67~ // - 課題なのですが検定統計量の使い方などがわかりません。 -- [[はえ]] &new{2022-07-28 (木) 23:30:52}; - 課題なのですが検定統計量の使い方などがわかりません。 -- [[はえ]] &new{2022-07-29 (金) 00:17:02}; #comment **剤投与状況の見える化のためのデータ加工 [#pb421e47] >[[shima]] (2022-06-02 (木) 10:57:46)~ ~ 患者毎の薬剤投薬状況の一覧を、見やすい形で提示したいのですが、データ加工方法がわからず悩んでいます。~ データの加工方法についてアドバイス頂けませんでしょうか？ネットや本Wikiで検索してもわからなくて・・・~ ~ 元のデータ~ 患者ID, 薬剤名, 投与量, 診断日, 投与日(診断日からの経過日数)~ 1111, drugA, 3, 2022/1/1, day2~ 1111, drugA, 2, 2022/1/1, day4~ 1111, drugA, 5, 2022/1/1, day8~ 1111, drugB, 1, 2022/1/1, day1~ 1111, drugB, 2, 2022/1/1, day4~ 1111, drugB, 1, 2022/1/1, day5~ 1111, drugB, 3, 2022/1/1, day6~ 1112, drugA, 7, 2022/1/7, day5~ 1112, drugA, 5, 2022/1/7, day8~ 1112, drugA, 2, 2022/1/7, day32~ ・~ ・~ ・~ ~ 最終的に提示したい形~ 患者ID, 薬剤名, 診断日, （投与日毎の投与量）day0,day1,day2,day3,day4,day5,day6,day7,day8,day9,day10,day11,day12,day13・・・~ 1111, drugA, 2022/1/1, 0,0,3,0,2,0,0,0,5,0,0,0,0,・・・・~ 1111, drugB, 2022/1/1, 0,1,0,0,2,1,3,0,0,0,0,0,0,・・・・~ 1112, drugA, 2022/1/7, 0,0,0,0,0,7,0,0,5,0,0,0,0,・・・・~ ~ 上記の最終的に提示したい形は、エクセルで見る場合に見やすいかな、と考えているものです。上記以外にも見やすい形式等のアドバイスをいただければ幸いです。~ // - うーん、ご自身で例示されているかたちで、とりあえず「患者毎の薬剤投薬状況の一覧」という目的は果たせているように（シロウト目には）見えますが…。~ 「ネットや本Wikiで検索しても」という行動をされているからには、なにか現状の提示方法では不都合や不便な点があるということでしょうか？具体的にそういった点があるのであれば明示していただくと、「なら、こうすればこの情報が見やすくなる」といった案も出やすいかもしれません。~ 個人的には、上記の例がどの程度、実データに即した例なのか不明なことも、他の提示方法を考えるうえでネックかなと思います。お書きになった提示方法だと表がすごく疎になってみづらいように感じますが、実際のデータでここまで間欠的に投薬を行なうことって少なそうで（科による？）、だとすれば例示された表が疎になっていること自体は大きな問題点でないのかもしれません。また薬剤ごとの投薬量や単位もバラバラなはずなので、整数値で並べて書ける（書くべき）ようなものなのかわかりませんし…。 -- &new{2022-06-02 (木) 13:48:14}; - コメントありがとうございます。患者数も薬剤名もかなり多かったので、何とか一覧で時系列で見やすくするために横長にしたかったというのがあります。含量と単位は薬剤名の中に入っており、まずは上記のような形式に変換する方法がわかれば、あとは元データからの加工方法を工夫すれば何とかなるかなと思っておりました（実際の元データは処方量と初報日数があり、日毎に処方量が提示され、それが処方日数分続くという表示が理想ですが、これも元データの加工次第かなと思っております）。 書き込んでから数時間また試行錯誤したのですが、pivot_widerで上記のような加工は可能そうでしたので、ひとまずは自己解決しました、すみません。 質問する際はもっと具体的に記載した方が良かったですね、お手数をお掛けいたしました、ありがとうございました。 -- [[shima]] &new{2022-06-02 (木) 15:14:17}; #comment **R-Tipsウェブサイト [#x9b8bbed] >[[山本]] (2022-05-18 (水) 12:10:24)~ ~ R-Tipsのウェブサイト（http://cse.naro.affrc.go.jp/takezawa/r-tips/r.html）が消えてしまったようですが、もしどこかに引っ越したのであれば、新しいURLをご存知の方、いらっしゃいませんか？~ // #comment **RjpWiki内で検索すると文字化けをする件 [#w56bc7ac] >[[かず]] (2022-04-03 (日) 08:56:09)~ ~ Rjpwiki内で特定の語句を検索しても文字化けで見れません(？ばかりが並びます)~ 対処方法をご教授いただけると助かります~ Microsoft EdgeやSafariいずれもです。~ // - 「検索しても」というのは、サイト上部のヘッダの「検索」機能のことですよね？わたしも試してみたところ、たしかに非Unicode文字みたいのが返ってきてるカンジでした。（ヘッダ部分は化けないので不思議です。）~ 取り急ぎ、Googleのサイト内検索は機能するので、そちらで検索していただいたらいかがでしょう？~ （Google検索時に「 site:http://www.okadajp.org/RWiki/」を付ける。）~ それとも、そうした「ただのサイト内検索」ではなく、あくまでPukiWiki機能での検索が必要な状況ですか？ -- &new{2022-04-03 (日) 09:10:46}; - ありがとうございました！！試してみます！ -- [[かず]] &new{2022-04-03 (日) 10:34:09}; #comment **openxlsx パッケージによる Excelファイルの編集 [#fc5dac4b] >[[ずっと初心者]] (2021-12-28 (火) 11:36:13)~ ~ Windows10 で Rを使っています。~ openxlsx パッケージを使って Excelファイル の読み書きをしたところ，~ 修正後のファイルでは 修正前に入力してあった日本語の文字列が化けてしまいました。~ 対処法をご存知の方がいらっしゃいましたら ご教示ください。~ options(encoding="UTF-8") wb <- openxlsx::loadWorkbook("sample1.xlsx") openxlsx::addWorksheet(wb, "A new worksheet") openxlsx::writeData(wb = wb, sheet = 1, x = "え", startCol = 8, startRow = 1) openxlsx::saveWorkbook(wb, "sample1_mod.xlsx", overwrite = TRUE) &ref(sample1.xlsx); &ref(sample1_mod.xlsx); >使用環境は下のとおりです。~ > sessionInfo() R version 4.1.2 (2021-11-01) Platform: x86_64-w64-mingw32/x64 (64-bit) Running under: Windows 10 x64 (build 19043) Matrix products: default locale: [1] LC_COLLATE=Japanese_Japan.932 LC_CTYPE=Japanese_Japan.932 [3] LC_MONETARY=Japanese_Japan.932 LC_NUMERIC=C [5] LC_TIME=Japanese_Japan.932 attached base packages: [1] stats graphics grDevices utils datasets methods base loaded via a namespace (and not attached): [1] compiler_4.1.2 tools_4.1.2 Rcpp_1.0.7 stringi_1.7.6 [5] zip_2.2.0 openxlsx_4.2.5 // - 直接の解決になっていなくて恐縮なのですが、当方のUbuntu16.04 ＆ R ver4.1.2環境ですと、上記ご提示のコードを実行しても文字化けは生じませんでした。（シートが増え、「Sheetあ」のH1のセルに「え」が書き込まれたxlsxファイルが保存される。）とりあえずopenxlsxパッケージそのものが本質的におかしいというワケでは無いようです。~ 関係あるかは定かでないですが、当方の環境だとlocaleまわりがUTF-8になっています。質問者さんの環境では932系になっているようですので、そのあたりの環境変数をいじってみたらいかがでしょう？ > strsplit(Sys.getlocale(), ";")[[1]] [1] "LC_CTYPE=ja_JP.UTF-8" "LC_NUMERIC=C" [3] "LC_TIME=ja_JP.UTF-8" "LC_COLLATE=ja_JP.UTF-8" [5] "LC_MONETARY=ja_JP.UTF-8" "LC_MESSAGES=ja_JP.UTF-8" [7] "LC_PAPER=ja_JP.UTF-8" "LC_NAME=C" [9] "LC_ADDRESS=C" "LC_TELEPHONE=C" [11] "LC_MEASUREMENT=ja_JP.UTF-8" "LC_IDENTIFICATION=C" 取り急ぎの参考情報まで。 -- &new{2021-12-28 (火) 19:58:33}; - ご返信ありがとうございます。 Sys.setlocale("LC_ALL", "ja_JP.UTF-8") などを試してみましたが，「ロケールを "ja_JP.UTF-8" に設定せよとの OS のレポート要求は受け入れられません 」という警告メッセージが出てしまい，途方に暮れております。 -- [[ずっと初心者]] &new{2022-01-04 (火) 11:20:58}; - なるほど...。そういえば以前Windows版ユーザだったころ、Rでlocaleまわりをいじろうとして意外と苦戦したことを思い出しました。（その際は日付の計算や曜日表記に関するプログラムで、英語に統一してどうにかしました。）ちょっとこちらでもWindows環境を探して試してみようかと思います。 -- &new{2022-01-04 (火) 12:15:21}; - 上記と同じ者です。たしかにWindowsだと、質問者さんのコードを実行しても失敗しますね。当方の環境だと、文字化けしたxlsxファイルができるのではなく、そもそも実行時にエラーとなりましたが...。不正なマルチバイト文字がみつかったときに、エラーで実行を止めるか、無理やり突き進むかの違いかなにかでしょう。 > options(encoding="UTF-8") > wb <- openxlsx::loadWorkbook("sample1.xlsx") 警告メッセージ: readLines(x, warn = FALSE, encoding = "UTF-8") で: 入力コネクション 'C:/Users/xxx/AppData/Local/Temp/RtmpXXXXXX/fileXXXXXXXXXXXX/xl/theme/theme1.xml' に不正な入力がありました ちなみにoptionsでUTF-8を指定しなかった場合でも、エラーは違いますが本質的には同様の理由でコケる感じですね。 -- &new{2022-01-04 (火) 12:22:56}; > wb <- openxlsx::loadWorkbook("sample1.xlsx") substring(x[ind], so, eo) でエラー: '<83><e3><82>ｯ"/>' に不正なマルチバイト文字があります - 理屈はイマイチ納得がいきませんが、とりあえず対処法としてlocaleをEnglishにしてみたら解決しました。（大いに不本意な感じですがw） > Sys.setlocale("LC_ALL", "English") [1] "LC_COLLATE=English_United States.1252;LC_CTYPE=English_United States.1252;LC_MONETARY=English_United States.1252;LC_NUMERIC=C;LC_TIME=English_United States.1252" > strsplit(Sys.getlocale(), ";")[[1]] [1] "LC_COLLATE=English_United States.1252" [2] "LC_CTYPE=English_United States.1252" [3] "LC_MONETARY=English_United States.1252" [4] "LC_NUMERIC=C" [5] "LC_TIME=English_United States.1252" > wb <- openxlsx::loadWorkbook("sample1.xlsx") > openxlsx::addWorksheet(wb, "A new worksheet") > openxlsx::writeData(wb = wb, sheet = 1, x = "え", startCol = 8, startRow = 1) > openxlsx::saveWorkbook(wb, "sample1_mod.xlsx", overwrite = TRUE) 結果、Linuxで出力されたのと同様の、文字化けのないsample1_mod.xlsxファイルができました。ただし注意として、その場合でもoptionsでUTF-8を指定すると失敗するようになります。（よくわかりませんw）~ 結局は日本語の932系エンコーディングが悪いことをしているっぽいのですが、localeの設定はOS依存なので、Linuxの設定をそのままWindowsでしようとしても弾かれたのでしょう。なので、とりあえあず無難な英語にSys.setlocaleして、あとは余計なoptionsとかをせずにopenxlsxを呼び出せば動いてくれる。。。ということで（納得はいきませんが）取り急ぎはどうでしょうか？ -- &new{2022-01-04 (火) 12:32:50}; - 色々と試してくださり ありがとうございます。こちらの環境でも optionsでUTF-8を指定して ['theme1.xml' に不正な入力が...」という警告が出ます。theme1.xml の中身を見ると，フォント名などに日本語やハングル文字が使われているようです。警告が出るだけで実行はできているようなんですが，出来上がる xlsxファイルはフォント名などを含めて日本語が文字化けしてしまうので困っています。optionsでUTF-8を指定しないと，「不正なマルチバイト文字があります」というエラーが出て実行ができません。 -- [[ずっと初心者]] &new{2022-01-04 (火) 12:46:33}; - 解決法をご提示くださり ありがとうございます。試してみます。 -- [[ずっと初心者]] &new{2022-01-04 (火) 12:47:47}; - Sys.setlocale("LC_ALL", "English") で うまくいきました！ありがとうございます！！！ -- [[ずっと初心者]] &new{2022-01-04 (火) 12:52:51}; - そちらでもうまく動いたようで良かったです。安心しました。わたしもopenxlsxという、いつか役に立ちそうなパッケージを知る機会になり勉強になりました。（私自身はExcelは使いませんが、共同研究者でどうしても使ってくる向きもみえられますので。） -- &new{2022-01-04 (火) 12:59:15}; - optionsでUTF-8を指定すると失敗するというご指摘の通り，その後，以下のようなコマンドを UTF-8 で書いて test_openxlsx.r として保存し，source("test_openxlsx.r", encoding="utf-8") を実行したら失敗しました。同様のコマンドで文字コードが sjis の "test_openxlsx_sjis.r" を作成して source("test_openxlsx_sjis.r") を実行したら うまくできました。-- [[ずっと初心者]] &new{2022-01-04 (火) 13:09:32}; Sys.setlocale("LC_ALL", "English") wb <- openxlsx::loadWorkbook("sample1.xlsx") openxlsx::addWorksheet(wb, "A new worksheet") openxlsx::writeData(wb = wb, sheet = 1, x = "え", startCol = 8, startRow = 1) openxlsx::saveWorkbook(wb, "sample1_mod.xlsx", overwrite = TRUE) - ぬーん、ややこしいもんですねぇ（笑）勝手なイメージとして、言語設定や文字エンコーディングまわりには「OSがユーザ向けフロントエンドとしてもってる設定」「OSのターミナル（コマンドプロンプト）とかが根っこの部分でもってる設定」「Rのもってる設定」「ファイル自体の属性」が絡んでいて、そのせいで挙動の理解が難しいのかも、と思っています。ファイル入出力やコマンドの受け渡しが起こる際に、それがどのレベルからどのレベルに対して起こっていて、その際に参照されるencodingやlocaleはどのレベルのなのか、といったような...。Windowsですら古いエンコーディングをどうにかする方向に動いている（のか？）っぽい昨今、OSの根幹レベルでUnicodeネイティブな時代がくれば、もっとスッキリ理解しやすくなる日が来るのかもしれませんね。 -- &new{2022-01-04 (火) 13:26:45}; #comment **反イメージ相関係数について [#c4d54729] >[[ともちゃん]] (2021-12-04 (土) 15:43:31)~ ~ Rで反イメージ相関係数を計算したいのですがわかりません。~ SPSSでの反イメージ相関係数の出力結果のような感じなのをRでもできるのでしょうか？ご教示いただけましたらありがたいです。~ // - 不勉強で反イメージ相関係数という言葉をはじめてお聞きしました。一部の統計処理ソフトのみに特有の用語のように見受けられますが、何のために、どのように算出される指標なのでしょうか？背景となっている数学的な理屈やその出典を示していただければ、それをプログラム的に実装してくれる回答者も出やすいかもしれません。逆に、数学的にどうやって算出されるものなのかがブラックボックスだったら、「Rで」に限らずどのプログラミング言語でも実装するのは困難であるように思います。 -- &new{2021-12-06 (月) 08:58:51}; - http://aoki2.si.gunma-u.ac.jp/lecture/PFA/pfa8.html 偏相関係数の符号を逆転したものは「反イメージ相関係数」と呼ばれ，この値がゼロに近いときは因子分析が有効であることを示すが，そうでない場合には得られたデータに対して因子分析を適用するのは不適切であることを意味する。 -- [[ ]] &new{2022-04-12 (火) 10:57:05}; - でしたらそのまま、偏相関係数を計算してマイナス符号をつけられたらよいのではないでしょうか？「SPSSでの反イメージ相関係数の出力結果のような感じなの」と仰々しくおっしゃられたので、何かSPSS特有の計算出力なのかと感じられましたが、そうではないのですねぇ。 -- &new{2022-04-13 (水) 07:26:51}; #comment **パッケージのインストールができない [#g0fbfb35] >[[もいくんとぅ]] (2021-08-21 (土) 21:38:06)~ ~ Rコマンダー、EZR、その他の必要なパッケージがインストールできないので、インストールできる方法をお教えして頂ければ有り難いです。~ ~ 使用OSは、Mac Book Pro 15inch Mid 2009 Mac 、Mac OS X Lion 10.7.5~ Rのバージョンは、R-3.2.1-snowleopard.pkg~ ~ ①CRANのホームページからインストーラーをダウンロードしインストール済みです。~ ②XQuartzプロジェクトからX11システムをインストール済みです。 ~ ③Rというアプリケーションを開き、install.packages("RcmdrPlugin.EZR", dependencies=TRUE)と入力、returnキーを押すが、下記エラーメッセージが表示されてしまいます。~ ~ 【警告: リポジトリー https://ism//cran.ism.ac.jp/src/contrib に対する牽引にアクセスできません。警告: リポジトリーhttps://ism//cran.ism.jp/bin/ Macosx/contrib/3.2に対する牽引にアクセスできません。警告メッセージ：package ’Rcmdrplugin.EZR'is Not available(for R version3.2.1】~ ~ CRANミラーをJapan(Tokyo),Japan(Yonezawa)に変更して試してみましたが、上記同様のメッセージが表示されてしまします。~ ~ 解決方法をお願い致します。~ // #comment **datasetsがありません [#g8a6f17a] >[[Rビギナー]] (2021-08-20 (金) 01:35:36)~ ~ コマンダーにおいて、パッケージ内のデーター→アタッチされたパッケージからデータセットを読み込む を行っても、CarDataとsandwichしかなく、datasetsが表示されません。どうしたらいいでしょうか、お教え願います。~ // - 初歩的な問題でした。テキスト付属のCDROMからアクセスできました。 -- [[Rビギナー]] &new{2021-09-16 (木) 13:32:39}; - そういう問題でしたか。解決されてよかったです。ご自身で解決報告をしていただけるのも、閲覧者目線としてはありがたいです。普段、計算関連の投稿ではできるだけ早く回答をつけたりしているのですが、RコマンダーやGUIは使わないもので、こういった系統のご質問にはご協力する術がなくROM専と化していました。 -- &new{2021-09-16 (木) 19:26:06}; - コメントありがとうございました。後からわかったのですがRコンソールでlibrary(Rcmdr)をロードさせると、datasetsが反映されるとわかりました。 -- [[Rビギナー]] &new{2021-09-18 (土) 14:18:55}; #comment **多変量正規分布の積分について [#ide39545] >[[よしの]] (2021-07-07 (水) 15:59:48)~ ~ 多変量正規分布の積分を行いたいのですが，使用できる関数がわかりません．~ ~ 正規分布（1変量）の場合，~ f_n <- function(x) dnorm(x, myu, sigma) integrate(f_n, -Inf, Inf) で計算でき，このような例はいくつかのサイトに示されていました．~ この1変量の場合と同様に，多変量でも積分を行いたいです．~ ~ 実際に試したいのは，3変量の場合で~ f_mv <- function(x1, x2, x3) dmvnorm(cbind(x1, x2, x3), mean=Mu1, sigma=Si1)~ と，関数を定義して~ これをx1, x2, x3について積分したいのですが，integrateを用いた場合の設定がわかりません．~ ~ 他の方法でも良いので，なにか情報いただければと思います．~ よろしくお願いします．~ ~ RStudio ver 1.3.1093~ // - cubatureパッケージを使うのが簡単そうですが、数値積分の用途でRを使ったことがないので自信はないですね……。 library(cubature) library(mvtnorm) Mu <- c(0, 1, 2) Sigma <- diag(length(Mu)) f_mv <- function(x) { dmvnorm(x, mean = Mu, sigma = Sigma) } hcubature(f_mv, c(-3, -2, -1), c(1, 2, 3)) パッケージのヘルプとVignetteを読んで使ってみてください。（わたしは斜め読みしただけなので。）どの程度の精度が必要かにもよるかもしれません。hcubatureはC言語で書かれた関数実体へのラッパーとのことなので、必要があればそちらのコードを確認するべきかと思います。 -- [[しろうと]] &new{2021-07-08 (木) 08:12:13}; - 情報ありがとうございます．いくつかの分布に対して処理を試してみました．多変量関分布のパラメタをデータから推定したものが，きちんと積分して１になるかをチェックしたいというのが目的でしたので，その用途としては十分です．ありがとうございました． -- [[よしの]] &new{2021-07-09 (金) 12:34:07}; #comment **プロビット分析の仮定 [#sf6abc04] >[[まるせん]] (2021-06-29 (火) 18:13:37)~ ~ 被説明変数に0,1の指標、説明変数に４項目程度の連続データでプロビット分析を行っています。データ数は400程度となります。~ Rでglmを使って無事モデルは構築できたのですが、右辺の被説明変数（標準正規分布の引数の部分）が標準正規分布に従っていません。（分布の形も正規分布ではなく、平均、分散もそれぞれ０，１ではありません。）~ そもそも右辺の被説明変数には標準正規分布の仮定を置いていないと思っていますが、私の考えがあっているか教えてください。前提は誤差項の正規性ですよね。~ これ読んどけ、というレポートを教えていただけるだけでもうれしいです。~ ~ お願いします。~ // - 右辺の被説明変数（？）云々のあたりがよく分かりませんが、プロビット分析の理屈自体としては「連続的な変数（群）からカテゴリカルな反応を説明するとき、両者をつなぐS字状の関数関係として正規分布の累積分布を使ってる」というだけではないですか？ -- &new{2021-06-30 (水) 07:40:30}; #comment **両対数目盛での回帰直線 [#xe0feb6d] >[[つぐ太郎]] (2021-06-19 (土) 15:47:10)~ ~ 年齢と死亡率が普通メモリグラフでは曲線ですが、両対数目盛グラフでは直線になることを示したいのですが、その時の回帰直線がうまく書けません。試したプログラムは下記です。~ age <- c(20, 25, 30, 35, 40, 45, 50, 55, 60, 65, 70, 75, 80, 85)~ mortality <- c(2, 5, 10, 20, 25, 37, 60, 80, 133, 190, 280, 410, 420, 380)~ plot(mortality ~ age, log = "xy")~ If <- glm(mortality ~ 1 + age)~ cf <- If$coefficients~ curve(cf[1] + cf[2] * x, from = 20, to = 80, col = "blue", cex=1.5, lwd=2, add = T)~ // - まぁ、出来合いの便利な関数や機能でやる方法もあるでしょうが。。。 -- &new{2021-06-20 (日) 13:37:58};~ 「両対数グラフで直線になる」ということは「ageのlogとmortalityのlogが線形関係にある」ということですよね？なので、 rslt <- lm(log10(mortality) ~ log10(age)) summary(rslt) としてやればよいわけです。 cf <- coef(rslt) par(mfrow = c(1, 2)) for (l in c("", "xy")) { plot(mortality ~ age, log = l) lines(age, 10 ^ (cf[1] + cf[2] * log10(age))) } - 明快にありがとうございます -- [[つぐ太郎]] &new{2021-06-22 (火) 23:36:34}; #comment **Rtoolsのパスが開けません [#x9a472fd] >[[KF]] (2021-04-10 (土) 20:44:50)~ ~ 以下を参考にrstanをダウンロードしようとしています。~ http://www.eeso.ges.kyoto-u.ac.jp/emm/materials/bayesian/stan_install~ rstan自体はダウンロード出来たのですが、rstanを使用するためにインストールしたrtoolsのパスを開こうとした際以下のエラーが出ました。~ ~ コード~ writeLines('PATH="${RTOOLS40_HOME}\\mingw64\\bin;${RTOOLS40_HOME}\\usr\\bin;${PATH}"', con = "~/.Renviron")~ 警告メッセージ: ~ file(con, "w") でエラー: コネクションを開くことができません 追加情報: 警告メッセージ: file(con, "w") で: ファイル 'C:/Users/*****/OneDrive/??????/.Renviron' を開くことができません: Invalid argument (*****＝ユーザー名) ~ 上のリンクでは「RでHomeとして使っているディレクトリ（デフォでgetwd()で表示されるパス）に，.Renviron というファイルが出来ているはず．」と書かれていましたが、そのようなファイルは作成されていないようでした。そのためエラーメッセージが発生したとみられます。~ ~ Rtoolsの設定がうまくいっていないせいか、rstanを起動しても以下のエラーメッセージが出ます。~ プロシージャエントリポイント EXTPTR_PTRがダイナミックリンクライブラリC:\Users\*****\R\win-library\4.0\rstan\libs\x64\rstan.dllから見つかりませんでした（*****=ユーザー名）~ ~ 操作環境は以下の通りです。~ R version 4.0.0 (2020-04-24)~ Platform: x86_64-w64-mingw32/x64 (64-bit)~ Running under: Windows 10 x64 (build 19042)~ ~ 解決方法をご教示いただけると幸いです。~ よろしくお願いいたします。~ // - 直接の解決策でなくすみません。~ 知り合いのPCをセットアップした際、「起動ディレクトリがOneDriveの配下だと.Rprofileファイルが読み込まれない」という問題にでくわしました。OneDriveまわりのWindowsの仕様は、R関連では邪魔をすることが色々あるようです。（検索してみてください。）管理者権限で起動していれば問題にならない場合もあるとか。~ .Rprofile問題の際は、OneDriveでないユーザホーム直下にファイル実体を置く、みたいな方法で対処した気がしますが、あまり定かではありません。OneDriveが無効な場所で.Renvironを置ける場所は他にないですか？あるいは古典的で泥臭いですが【システムのプロパティ】→【環境変数】でPATHをR外で設定しておくとか。 -- &new{2021-04-11 (日) 14:17:42}; - コメント頂きありがとうございます。OneDriveではない場所に新規のプロジェクトを作成しました。writeLines('PATH="${RTOOLS40_HOME}\\mingw64\\bin;${RTOOLS40_HOME}\\usr\\bin;${PATH}"', con = "~/.Renviron")と入力したところ、前回と全く同じエラーが出てしまいました。file(con, "w") でエラー: コネクションを開くことができません 追加情報: 警告メッセージ: file(con, "w") で: ファイル 'C:/Users/*****/OneDrive/??????/.Renviron' を開くことができません: Invalid argument (*****＝ユーザー名) 何が問題なのか全く分からず途方に暮れております。どうぞよろしくお願い致します。 -- [[KF]] &new{2021-04-12 (月) 18:21:14}; - R, Rstudio, Rtoolすべてアンインストールし、再インストールしたらrstanが正常に使えるようになりました。保存場所をOne Driveではないフォルダにしたことが効いたのかもしれません。ありがとうございました。 -- [[KF]] &new{2021-04-12 (月) 20:37:54}; #comment **read.table()で部分的に文字化けしてしまう [#jc87185b] >[[Gorou]] (2021-04-07 (水) 20:41:39)~ ~ 質問です。"UTF-8"ファイルの読み込み時に文字化けしてしまい困っています。~ （Wikiの記法に不慣れな点はお許し願います。)~ ~ read.table()での文字化けは非常にありがちな質問ですが、~ ~ "read.table()"で読み込んだデータが"「部分的に」文字化け"してしまい、~ ~ 困っています。~ ~ 具体的には次のような症状です。~ 以下のデータファイルを読み込みたいのですが、~ [data.tsv] (TSV形式、UTF-8、BOMなし） Age Group Class DBH_mean DBH_sd 31 すべて 全立木 17.9 4.34 31 すべて 生立木 17.9 4.34 31 すべて 残存木 17.9 4.34 36 すべて 全立木 20.0 4.92 36 すべて 生立木 20.0 4.92 36 すべて 残存木 20.0 4.92 41 すべて 全立木 21.9 5.59 41 すべて 生立木 21.9 5.56 41 すべて 残存木 21.9 5.56 ~ >df <- read.table("test.tsv",sep="\t",header=TRUE,skip=0, fileEncoding="UTF-8") >head(df) ~ とすると、~ ~ Age Group Class DBH_mean DBH_sd X 1 31 すべて \u5168\u7acb\u6728 17.9 4.34 NA 2 31 すべて \u751f\u7acb\u6728 17.9 4.34 NA 3 31 すべて \u6b8b\u5b58\u6728 17.9 4.34 NA 4 36 すべて \u5168\u7acb\u6728 20.0 4.92 NA 5 36 すべて \u751f\u7acb\u6728 20.0 4.92 NA 6 36 すべて \u6b8b\u5b58\u6728 20.0 4.92 NA ("\"はバックスラッシュ） ~ となり、「Group」のカラムは文字化けしないのに、「Class」のカラムは文字化けしてしまいます。完全に文字化けしているのならば文字コードの間違いだと納得できるのですが、~ このような中途半端な文字化けは、いろいろググって見たのですがさっぱり分かりません。~ ~ また、DBH_sd の右に "X"のカラムができてしまうのも不可解です。DBH_sdのカラム、4.34, 4.34, 4.34, ... は数値の直後で改行されています（つまり元のファイルにはヘッダだけのカラムは無い）。~ ~ TSVをCSVにしたり、ファイル形式をShift-JISにして fileEncoding='Shift_JIS'としてみても、結果はほぼ同じで途方に暮れております。~ おまじないにと,stringAsFactorsをTRUEにしたりFALSEにしても変わりません(変わるとも思えませんが）。~ ~ 環境はWindows10 64bit Pro、Rのバージョンは64bitの4.0.4、Rstudioは1.4.1106。症状はRGuiでもRStudioでも同じです。~ ~ どなたか解決方法をご教示いただけましたら幸いです。何卒よろしくお願いします。~ // - とりあえず列が増えるのは、単に行末の余計なタブのせいでしょう。質問者さんが貼られたデータ例にも、行末にタブが入っています。タブのうしろに何もデータがないので、read.table時に自動で列がつくられたものの、値はすべてNAになっています。 ~ ~ つぎにClassの列の現象ですが……文字化けというよりは、UTF-8のコードに展開されているというカンジでしょうか。数値をUTF-8として読み取ればもとの「全立木」などの文字列が読み取れますので、データとして「化けている」というのともちょっと違う気がします。~ ~ なぜそんなことになってるかについては……難しいですねぇ。可能性として、以下のようないろんなパターンが考えられるかと思います。質問者さんが行末のタブに気づいていらっしゃらなかった（数値後即改行だと思われていた）件もあわせると、最後のような「じつはR外の問題」みたいな可能性も無きにしも非ずではないかな、と。（僭越ですが。）その現象が起こるminimal exampleのファイルをご提示いただければ、どこに問題があるかの切り分けが多少楽になるかと思うのですが。 --その列のデータだけ何か理由があって、read.tableがおかしな挙動をした。 --read.tableはできてるんだけど、それを表示するメソッドprint.data.frameがおかしな表示の仕方をした。（たとえばその列だけ特殊なclass属性がついちゃってて、print時にそれが反映されてしまった。） --じつはもとのデータもUTF-8のコード直書きみたいな変な形式になってるのだけど、（R外でそれを表示して）確認に使ってるソフトが自動的に文字列として表示してくれてるので気づけていない。Rはそれを愚直に読み込んでるだけ。 -- &new{2021-04-08 (木) 13:45:15}; - exampleファイルです -- [[Gorou]] &new{2021-04-09 (金) 13:09:14}; - 済みません、操作ミスです（削除方法を勉強中です…） -- [[Gorou]] &new{2021-04-09 (金) 13:16:44}; -- 列が増えるのは、余計な半角スペースが入っていたためでした。削除したら解決しました。 -- [[Gorou]] &new{2021-04-09 (金) 13:19:13}; -- オリジナルのデータは、RStudio、秀丸エディタ、notepad++のいずれでもUTF-8として認識、表示されています。 -- [[Gorou]] &new{2021-04-09 (金) 13:22:20}; -- (ubuntuで）nkf -g test.tsv とするとUTF-8が帰ってきますし、"cat test.tsv" でも正しく表示されます。（minimal exampleのファイルの提示方法が分かりません。箇条書きですみませんが、よろしくお願いします。） -- [[Gorou]] &new{2021-04-09 (金) 13:27:07}; -- 度々済みません。根本的な解決策ではないのですが、"read.table"の代わりに"tidyverse"を使ったら文字化けせずに読み込めました。しかし、以前作ったデータとスクリプトの組み合わせが使えないままです。引き続き、情報等ありましたらよろしくお願いします。(Wikiの記法、勉強します…。） -- [[Gorou]] &new{2021-04-09 (金) 14:14:50}; - 確認ですが、このページに添付されたtest_utf8.tsvというファイルは、質問者さんの環境では前記のような問題の再現するファイルなのですか？わたしの環境だと、その「文字化け（と呼んでいらっしゃる現象）」なしでちゃんと読み込めるようです。 > df <- read.table("test_utf8.tsv", sep = "\t", header = TRUE, skip = 0, fileEncoding = "UTF-8") > head(df) Age Group Class D_mean H_mean N G V 1 31 すべて 1)間伐前 17.9 12.3 1535 40.95 260.0 2 31 すべて 2)間伐後 17.9 12.3 1535 40.95 260.0 3 36 すべて 1)間伐前 20.0 14.9 1535 51.33 393.1 4 36 すべて 2)間伐後 20.0 14.9 1535 51.33 393.1 5 41 すべて 1)間伐前 21.9 16.6 1535 61.40 517.8 6 41 すべて 2)間伐後 21.9 16.6 1530 61.36 517.6 なおこちらは、R version 4.0.4 on Ubuntu 16.04 LTSをbash上で起動しています。 -- &new{2021-04-09 (金) 17:09:52}; - （最初に提示されたデータと列などが違うので、別のデータなのかな……もしunpublishedで外に出たらまずいデータとかだったら、消すなり差し替えるなりなさってくださいね。） -- &new{2021-04-09 (金) 17:18:14}; - 確認ですが、このページに添付されたtest_utf8.tsvというファイルは、質問者さんの環境では前記のような問題の再現するファイルなのですか？ > はい、その通りです。データを確認していただき、ありがとうございます。ちょっと書き換えましたがほぼ同様のデータです。 理由はさっぱり分からないのですが、Rを4.0.4から4.0.5にアップデートしたら文字化けしなくなりました(RStudioもRGuiも）。 違う点としては4.0.5は管理者として、4.0.4は一般ユーザとしてインストールしたことです。 （管理者だとC:\Program Files\R\R-4.0.5に、一般ユーザだとC:\Users\hogehoge\Documents\R\R-4.0.4にインストールされていました。） それにしても、元のデータは > すべて 全立木 であるのに「文字化け」状態だと > すべて \u5168\u7acb\u6728 と部分的に化けてしまうと言うのはさっぱり訳が分かりません。 いずれにしても、これまでのデータやスクリプトが問題なく使えるので安心しました。 助言をくださったみなさん、どうもありがとうございます。 （ちなみにtidyverseのread_csvも、再確認したところ4.0.5にアップデートする前にはなぜか文字化けしてました。こっちも何が何やら。） - （返信の編集方法はこれで良かったのかな？みなさん、ありがとうございます。） -- [[Gorou]] &new{2021-04-12 (月) 16:55:33}; #comment **シミュレーションで確率を見積もりたい [#lba761df] >[[HO0217]] (2020-12-29 (火) 21:28:39)~ ~ 「赤札と青札がある。確率1/2でどちらかを引く。~ これを100回繰り返し、赤札を55回以上引く確率を求めよ。」~ という問題を解くため、正規分布に近似して確率を求めると、~ 15.87%でした。~ (これは私の計算結果で、答えは知りません。あってるとは思いますが、間違っているかもしれません)~ ~ ~ 近似ではなく、二項分布で計算してみようと思い、~ 以下のコードで確率を計算すると、約18.41%(出力は0.1841008)で、~ 近似とはやや違う結果となりました。~ tmp_nu<-0~ n_nu<-100~ k_nu<-55~ for(i in k_nu:n_nu){~ tmp_nu<-( choose(n_nu,i)*(1/2)^(i+(100-i)) )+tmp_nu }~ tmp_nu~ ~ ~ 次に、シミュレーションで確率を見積もるために以下のコードを実行しました。~ 出力は0.13554となりました。~ ~ simu_0_nu<-100 #100枚引く。~ simu_1_nu<-1000000 # 何セット試すか。~ result_0_ve<-rep(NA,simu_0_nu)~ result_1_ve<-rep(NA,simu_0_nu)~ result_2_ve<-rep(F,simu_1_nu)~ set.seed(123)~ for(i in 1:simu_1_nu){~ for(j in 1:simu_0_nu){ result_0_ve<-rnorm(simu_0_nu,0,1) #平均0の正規分布は確率1/2で0を越える数値を出す result_1_ve<-result_0_ve>0 #0を越えたら真 } result_2_ve[[i]]<-sum(result_1_ve)#真の数を数える }~ ~ sum(result_2_ve>55)~ sum(result_2_ve>55)/simu_1_nu~ ~ ~ ちなみに、set.seed(123)を~ set.seed(777)にすると0.135447~ set.seed(789)にすると0.136265~ set.seed(111)にすると0.136006~ set.seed(222)にすると0.135414~ となりました。~ どうやら偶然ではなく、私がシミュレーションに使ったコードでは、~ 他の計算よりも小さい値になっているようです。~ ~ 何が原因でシミュレーションと二項分布の計算でこんなにも違ってしまうのでしょうか？~ ~ わかる方がいらっしゃいましたら、教えてくださると助かります。~ よろしくお願いいたします。~ // - 55を含むかどうかの差ではないのですか？ -- &new{2020-12-30 (水) 08:27:58}; > set.seed(1) > N <- 1000000 > n <- replicate(N, {sum(rnorm(100) > 0)}) > sum(n > 55) / N [1] 0.135948 > 1 - pbinom(55, 100, 0.5) [1] 0.1356265 > sum(n >= 55) / N [1] 0.184649 > 1 - pbinom(54, 100, 0.5) [1] 0.1841008 - ありがとうございます。全く気が付きませんでした。 -- [[HO0217]] &new{2020-12-30 (水) 09:46:44}; #comment **所感 [#y895d2fe] >[[元Rユーザー]] (2020-12-27 (日) 22:33:22)~ ~ 全盛期とは比較にならないほどの過疎ぶりですね・・・時代は完全にpythonになったのだなと思います。~ // - ここは初心者に辛辣なコメントが多いのも一因かなぁと。slackのr-wakalangとかは賑わってますよ -- [[Rユーザー]] &new{2021-01-13 (水) 17:29:49}; - 「一時のような辛辣なコメントが減って落ち着いたので、良いんじゃないかなぁ？」と思っていたところです（笑） -- [[現Rユーザ]] &new{2021-01-15 (金) 09:18:03}; #comment **パネルデータ分析における警告メッセージ（重複エラー？） [#w167b8e1] >[[R初心者]] (2020-11-08 (日) 22:39:38)~ ~ パネルデータでplm分析を行ってますが、下記のメッセージが出て対応に苦慮しております。 ~ ~ 警告メッセージ: ~ 1: pdata.frame(data, index) で: ~ duplicate couples (id-time) in resulting pdata.frame to find out which, use e.g. table(index(your_pdataframe), useNA = "ifany") 2: is.pbalanced.default(index[[1]], index[[2]]) で: ~ duplicate couples (id-time) ~ 恐らく重複エラーなのかもしれませんが、下記を見ると年次などの下線部がエラーの様ですが（？）、削除はしたくなく対応方法を知りたいです。宜しくお願いいたします。~ ~ コード PBR 業種コントロール 年次 `1人当たり単体人件費`… 連結子会社数 売上高研究開発費 <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> 1 T1301 0.915 0.630 2019 5.87 26 0.0013 2 T1301 0.908 0.618 2018 5.88 24 0.00120 3 T1301 0.872 0.644 2017 6.18 24 0.0011 4 T1301 0.867 0.678 2016 6.15 25 0.00120 // #comment **急ぎで知りたいです。 [#k20d893a] >[[Rビギナー]] (2020-10-28 (水) 01:22:23)~ ~ SVM実装の際に、package(kernlab)を使用しています。ksvm関数を使用すると以下のようなエラーが出ます。これはどういう意味でしょうか？~ ~ SVM <- ksvm( data_Mic.dummy ~AveMaxTemp7+TP1yp+TN1yp+AveWind7, #モデル式~ data=data_Mic, kernel ="rbfdot", #カーネル関数の設定 kpar=list(sigma=0.2), #カーネル関数のパラメータ C=1, cross=10) .local(x, ...) でエラー: ~ No Support Vectors found. You may want to change your parameters 追加情報: 警告メッセージ: .local(x, ...) で: 強制変換により NA が生成されました // #comment **Rで複数の関数のセットをループしたい。 [#u1756c69] >[[R初心者]] (2020-10-20 (火) 12:55:09)~ ~ R初心者で色々わかっていないので、質問の仕方が間違ってたらすみません。~ Rで関数の処理を行いたいのですが、~ ~ 解析結果のデータは~ ~ dds <- 解析結果~ ~ に格納しました。~ その大量のデータの中で図示したい５つだけをリストアップし~ dataという下記のようなデータフレームを作りました。~ ~ name ID~ 1 aaaaa 001~ 2 bbbbb 002~ 3 ccccc 003~ 4 ddddd 004~ 5 eeeee 005~ ~ 処理したい関数が下記のように３つでpng→plotCount→dev.offでワンセットです。~ ~ 図示したいのが1つであれば~ png("aaaaa.png")~ plotCounts(dds, gene="001", intgroup="condition" )~ dev.off~ ~ でデータフレームの一行目のものが図として保存されます。~ ですが連続してdataの二行目、三行目...と自動でループ処理させたくて、~ forやlapplyなどを試したのですが、引数がnameとIDの二つということもあり、うまくプログラムを回せませんでした。~ ~ どなたかお詳しい方に上記の関数でpng→plotCount→dev.off→png→plotCount→dev.off→...を自動的に回す方法をご教授頂きたいです。~ ~ よろしくお願いします。~ // - データフレームの中身をループで取得してください for(i in 1:nrow(dds)){ name_ <- dds$name[i]; ID_ <- dds$ID[i] } ループの中でpng, plotCount, dev.offを順次呼び出せば動くと思います。 -- &new{2020-10-23 (金) 11:55:09}; - ありがとうございます。 -- [[R初心者]] &new{2020-11-06 (金) 17:23:08}; #comment **回帰直線 (RMA) の描き方と統計情報の出し方等について [#bec2344c] >[[うみねこ]] (2020-07-24 (金) 22:13:58)~ ~ R version 3.6.1 (2019-07-05)~ Platform: x86_64-w64-mingw32/x64 (64-bit)~ Running under: Windows 10 x64 (build 18363)~ ~ 回帰直線の描き方と、その統計情報の出し方について、教えて頂けないでしょうか。~ ９列×６行の表があり、その２列目（length）と６列目(weight)にlog変換をして、RMA法を使って回帰直線を描きたいです。gx.rmaというファンクションを見つけましたが、以下のコードではエラーが出て動きません。~ 出来たら、95％予測区間や95％信頼区間等についても出し方を知りたいです。~ どうぞご教授よろしくお願いいたします。~ ~ data<- read.table('body.txt', header=T)~ row.names(data)=data[,1]~ length <- setNames(log(data[,2]),data[,1])~ weight <- setNames(log(data[,6]),data[,1])~ model <-gx.rma(length,weight)~ model~ // - 実際の data がどういうものかわからないので，適当にいかのようなデータをつくってやってみたけど，エラーは出ませんでしたが。（質問の仕方がヘタな人がそろっちゃったなあ） > set.seed = 1234567 > data = as.data.frame(matrix(runif(54), 6, 9)) > data[,1] = letters[1:6] > row.names(data)=data[,1] > length <- setNames(log(data[,2]),data[,1]) > weight <- setNames(log(data[,6]),data[,1]) > model <-gx.rma(length,weight) Reduced Major Axis for length and weight Means = -1.05 -0.679 SDs = 0.7428 0.8641 Corr = 0.5295 N = 6 Fit = 28% SE 95% CLs Slope = 1.163 0.4028 0.1277 <-> 2.199 Intercept = 0.5422 0.518067 -0.7895 <-> 1.874 Accept hypothesis that RMA is (0,1) > model $n [1] 6 $mean [,1] [,2] [1,] -1.049833 -0.6789889 $sd [,1] [,2] [1,] 0.7428316 0.8640773 $corr [1] 0.5295494 $a0 [1] 0.5421989 $a1 [1] 1.163221 $sea0 [1] 0.5180675 $sea1 [1] 0.4028333 他の人が追試できるように質問すべし。 -- &new{2020-07-26 (日) 17:41:28}; - すみません。質問が下手で、大変失礼いたしました。データですが、一列目が、人の名前、一行目は身体の長さ、幅、重さ等になっています。再度実行してみましたら、以下のようなエラーが出てしまい、まだ出来ていません。ご教授よろしくお願いいたします。 gx.rma(length, weight) でエラー: 関数 "gx.rma" を見つけることができませんでした -- [[うみねこ]] &new{2020-07-26 (日) 22:19:06}; - > 関数 "gx.rma" を見つけることができませんでした&br;って，library(rgr) やってないんですか？あるいは，install.packages("rgr") さえもしていないとか。 -- &new{2020-07-27 (月) 07:34:13}; - install.packages("rgr") をしていないのが原因でした。確認したところ、gx.rmaのファンクションの記載があるドキュメントの最初のページにrgrと書いてありましたが、スルーしていたようです。無事に上述の結果が出ました！ありがとうございます。次は、この結果を（95％信頼区間と回帰直線）プロットしようとして、plot(model)と打ち込みましたが、xy.coords(x, y, xlabel, ylabel, log) でエラー: 'x' is a list, but does not have components 'x' and 'y'と出て失敗しました。さらに、以下のように打ち込んでみましたが、plot (length, weight,lengthlab=”Length(cm)”,weightlab=”Weight(g)”) エラー: 想定外の入力です in "plot (length, weight,lengthlab=・と出て失敗しました。最小二乗法と同じように作っているつもりですが、plotも出来ません。何度もすみませんが、再度ご教授どうぞよろしくお願いいたします。 -- [[うみねこ]] &new{2020-07-27 (月) 09:44:23}; - > 最小二乗法と同じように作っているつもりですが&br;どの最小二乗法？そもそも，gx.rma の plot method はないと思いますが。あなたが参照している（？）「gx.rmaのファンクションの記載があるドキュメント」に，どのように書いてありますか？書いてなくて，自分で適当に書いてもそれは無理です。lengthlab=”Length(cm)”, weightlab=”Weight(g)”) なんていう引数は plot() にはないのでは？&br;オンラインヘルプを見るとわかりますが，gr.rma の引数に ifplot というのがあるのに気づいていますか？ ifplot if a x-y plot of the independent variables is required set ifplot = TRUE. The plot is equi-scaled and the 1:1 line is added. If ifrma = TRUE ifplot will be set to TRUE. よくわからないまま，適当にプログラムを書いてもだめです。オンラインヘルプくらい読みましょう。それと，エラーメッセージを良く読めば，どこが間違えているか対処法もわかります（わかりやすく書かれているとはいえないでしょうが，それが問題解決の唯一の鍵ですからね） -- &new{2020-07-27 (月) 14:41:47}; - はい。gx.rmaには掲載されていません。if plotは、もちろん気づいていますし、英語の意味も取れますが、それの入力の仕方がわかりません。だから、代表的コマンドの、Plot(x,y)で何とかならないかと、ついでにabline(lm(y,x))も試してみました。おそらくこれらは、ｙ軸方向に残差を取るタイプの一番メジャーな最小二乗法だと思います。オンラインヘルプを読んでも、エラーメッセージを読んでも、読めるけれども、実力がなくて対処が出来ないというのが問題で、初心者打破のため、いろいろ作成して練習していますが、難しいです。 -- [[うみねこ]] &new{2020-07-27 (月) 16:23:18}; - > それの入力の仕方がわかりません&br;Usage の下に，関数の使い方（引数の種類と，その指定法）が書いてあるじゃないですか。 Usage gx.rma(xx1, xx2, x1lab = NULL, x2lab = NULL, log = FALSE, ifplot = FALSE, ifrma = FALSE, ifcoeffs = FALSE, ifform = FALSE, iftest = FALSE, ...) あなたの例の場合だと， gx.rma(length, weight, x1lab="Length(cm)", x2lab="Weight(g)", ifplot=TRUE, ifrma=TRUE) のようになる。 -- &new{2020-07-27 (月) 16:27:48}; - ありがとうございます。無事にplotされました。このUsageもですが、もしgx.rma(x,y,xlab=NULL, ylab=NULL...)であったら、もう少し私にはわかりやすいんです。使用例のコードを見ると、余計にUsageが解読し辛くなってしまう位のレベルですが、諦めずに頑張ろうと思います。ありがとうございます。 -- [[うみねこ]] &new{2020-07-27 (月) 17:47:53}; #comment **Treatment indicator ('Tr') must be a logical variable---i.e., TRUE (1) or FALSE (0) [#bdc9264a] >[[catdogcat]] (2020-07-24 (金) 18:20:51)~ ~ マッチングを試みましたが下記のエラーが出ました~ ~ Treatment indicator ('Tr') must be a logical variable---i.e., TRUE (1) or FALSE (0)~ ~ 比較する群の変数は因子としてあり 0.1 となっているのですが~ 解決策がわかりません。~ ご教授お願い致します。~ // - 実際にあなたが何をやったかが，全くわからない。なのに，解決策を授けられる人がいるだろうか？いや，いない。&br;「群の変数は因子としてあり 0.1」ということは，"must be a logical variable" という条件を満たしていないのではないかな？ -- &new{2020-07-26 (日) 17:46:21}; - 幼稚な質問であり申し訳ございません。確かにlogical variableになってなかったようでした。ありがとうございます。 -- [[catdogcat]] &new{2020-07-26 (日) 22:23:43}; #comment **Rcppで行列を多次元配列にfor文で代入する [#m2aaf0fb] > (2020-07-22 (水) 23:44:38)~ ~ 下記コードをpgm.cppに保存後、Rcpp::sourceCpp("pgm.cpp")でコンパイルし~ test(matrix(c(1,2,3,4),2,2))を実行したのですが1行1列の値が返ってこずエラーでRが落ちます。もし良ければご助言いただけないでしょうか。~ ~ ~ #include<Rcpp.h> #define MAX 100 using namespace std; // [[Rcpp::export]] int test(Rcpp::IntegerMatrix mtx){ int matrix[MAX][MAX]; for(int i=0; i < mtx.rows(); ++i){ for(int j=0; j < mtx.cols(); ++j){ matrix[i][j] = mtx(i,j); }} return matrix[0][0]; } // - やってみたけど，確かに1x1の行列は帰ってこないけど，R が堕ちるというようなことはなかったけど？ -- &new{2020-07-26 (日) 17:47:18}; #comment **くの字に曲がるデータプロットの特異点の算出？ [#b8b95798] >[[がんば]] (2020-07-22 (水) 08:35:12)~ ~ 昨年はおせわになりました。実験結果をグラフにすると直線状ではなく「くの字」状に曲がった結果となる場合もあるのですがこの時の特異点をＲで一発算出するような手段がないでしょうか？「データをソートして下限側から徐々に近似曲線を構築するのと上限側から構築するのをそれぞれ行っていって次の点がそれまでの近似曲線のバラツキ範囲外になる点を求める」ような操作をする（人の目視確認の場合はこんな手順だと思います）手段を意図してます。~ 動作環境：~ > sessionInfo() ~ R version 4.0.2 (2020-06-22)~ Platform: x86_64-w64-mingw32/x64 (64-bit)~ Running under: Windows 10 x64 (build 18363)~ Matrix products: default~ locale:~ [1] LC_COLLATE=Japanese_Japan.932 LC_CTYPE=Japanese_Japan.932 LC_MONETARY=Japanese_Japan.932~ [4] LC_NUMERIC=C LC_TIME=Japanese_Japan.932 ~ attached base packages:~ [1] stats graphics grDevices utils datasets methods base ~ loaded via a namespace (and not attached):~ [1] compiler_4.0.2 tools_4.0.2~ // - Googleで"bent line gression R"で検索すると、Rbentパッケージが出てきたけど、rbentfit()が吐くbpがもしかしたら、お望みのものでは -- &new{2020-07-22 (水) 22:13:43}; - http://aoki2.si.gunma-u.ac.jp/R/oresen.html 二本の直線による折線回帰 とかは？ -- &new{2020-07-26 (日) 17:17:55}; - ↑を元に検討してみます。ありがとうございます。 -- [[がんば]] &new{2020-07-29 (水) 08:41:50}; - RbentパッケージについてはhelpのExsamplesをいじってみて検討します。ありがとうごじます。 -- [[がんば]] &new{2020-07-29 (水) 09:01:05}; #comment **R studio でグラフ化した時に日本語文字化けが起こる [#bacaed2a] >[[りょう]] (2020-06-30 (火) 21:38:04)~ ~ Rによる多変量解析と言う本でRを動かしているのですが、グラフ化すると日本語の部分が□□□のようになります。~ よくわからないのが、"plot", "plotmeans"などでは日本語表示されるのですが、xyplot, histogramだと文字化けすることです。以下がその時のコードです。~ ~ >jhk <- read.csv('人事評価結果.csv')~ head(jhk, 3)~ |ID|性別|部署|年代|協調性|技能|知識|ストレス|~ |1|M|A部|中堅|70|65|71|53|~ |2|F|B部|熟練|45|62|70|46|~ ・・・・・・~ par(family = "HiraKakuProN-W3") #一応これ試してみました~ library(lattice, gplots)~ histogram(~ストレス, data = jhk, breaks = 20, type='count') #文字化けする~ xyplot(知識~技能|年代+部署, data=jhk) #文字化けする~ plotmeans(協調性~性別, data=jhk, p=0.95, ylim=c(49, 54)) #文字化けしない~ ~ ・macOS Catalina バージョン10.15.5~ ・R studio version 4.0.2 です。この前に使用していたバージョン3.6.3でも同様の状態でした。~ ・CSVファイルはUTF-8で保存されています。~ ・par(family = "HiraKakuProN-W3") なども試しました。~ ~ GoogleやCRAN, マニュアルなどで既に調べたのですがわからず、困っているのでご教授いただけると幸いです。~ // - ごめんなさい。コードの部分を青くしたり、表形式にしたつもりができていません汗 -- [[りょう]] &new{2020-06-30 (火) 21:39:55}; - お門違いの場所にアップされ，誰も返事をしない，もう見ていないかもしれないが一応書いておこう。&br;histogram, xyplot は lattice のものなので，日本語フォント指定が旨くいかない（なんか方法はあるはずだが，ネットをあさっても...）。&br;回避法としては&br;「画面に描かず PDF ファイルに保存する」とちゃんと描かれるっす。 pdf("test.pdf") # ひつようなら他のパラメータの設定も histogram(~ストレス, data = jhk, breaks = 20, type='count') #文字化しない xyplot(知識~技能|年代+部署, data=jhk) #文字化けしない plotmeans(協調性~性別, data=jhk, p=0.95, ylim=c(49, 54)) #文字化けしない dev.off() どんなもんだ？ -- &new{2020-07-14 (火) 14:04:24}; #comment **ggplot [#bc6d7d37] >[[ぼう]] (2020-05-14 (木) 12:25:32)~ ~ R4.00でggplotを使おうと思っていますが、エラーが出てしまいます。~ ~ 次のコマンドを実行すると、下のような結果が出てきてしまいます。 > install.packages("ggplot2", dependencies = TRUE) > library(ggplot2) > str(iris) 'data.frame': 150 obs. of 5 variables:~ $ Sepal.Length: num 5.1 4.9 4.7 4.6 5 5.4 4.6 5 4.4 4.9 ... $ Sepal.Width : num 3.5 3 3.2 3.1 3.6 3.9 3.4 3.4 2.9 3.1 ... $ Petal.Length: num 1.4 1.4 1.3 1.5 1.4 1.7 1.4 1.5 1.4 1.5 ... $ Petal.Width : num 0.2 0.2 0.2 0.2 0.2 0.4 0.3 0.2 0.2 0.1 ... $ Species : Factor w/ 3 levels "setosa","versicolor",..: 1 1 1 1 1 1 1 1 1 1 ... > ggplot(data=iris,aes(x=Sepal.Length,y=Sepal.Width))~ ggplot(data = iris, aes(x = Sepal.Length, y = Sepal.Width)) でエラー: 関数 "ggplot" を見つけることができませんでした // - ggplot2 がインストールされていない（OS はなんですか。Mac では R4.0 なら，コンパイルが必要なインストールが延々と実行されますが，それがちゃんと出来ているか。Windows ならそうでもないのかな）~ と，思ったが，Windows でもインストールはコンパイルが必要なようですね。で，それでも，ちゃんとインストールできましたが。~ あなたの場合，コンパイルに失敗していると言うことはありませんか。コンパイルが可能な環境が出来ていないのかな。普通は，そんなことないはずなのだけど。 -- &new{2020-05-15 (金) 21:38:52}; - ありがとうございます。Win10を使っています。library(ggplot2)でエラーが出ていないのに、ggplotを実行するとエラーが出てしまうので、原因が分かりません。。。 -- [[ぼう]] &new{2020-05-18 (月) 11:28:18}; - ありがとうございます。Win10を使っています。library(ggplot2)でエラーが出ていないのに、ggplotを実行するとエラーが出てしまうので、原因が分かりません。。。 -- [[ぼう]] &new{2020-05-18 (月) 12:09:38}; #comment **ccf [#y347335b] >[[show]] (2020-03-10 (火) 12:58:28)~ ~ ccf(rw,rw2)と入力した場合、どちらを基準にしてlagを取っているのでしょうか？~ // - ? ccf -- "The lag k value returned by ccf(x, y) estimates the correlation between x[t+k] and y[t]." -- &new{2020-03-10 (火) 17:55:37}; - rw2を基準としているということでしょうか。 -- [[show]] &new{2020-03-11 (水) 12:34:59}; - what? rw2? what do you mean?-- &new{2020-03-11 (水) 22:56:56}; - ccf(rw,rw2)の場合、rw２を基準にしてるということでしょうか？ -- [[show]] &new{2020-03-13 (金) 08:33:16}; - x[t+k] and y[t] -- &new{2020-03-13 (金) 10:12:16}; - なぜ自分で確かめもせず，不毛な質問を続けるのかなあ？~ 以下のようなテストデータを作ってみる。ベクトルの長さ n = 10 > set.seed(12345) > n <- 10 > x <- rnorm(n) > y <- rnorm(n) ラグ lag = 2，で ccf してみる > print(ccf(x, y, lag, plot=FALSE), digits=7) Autocorrelations of series ‘X’, by lag -2 -1 0 1 2 0.1416173 0.5327757 -0.2716881 0.2595251 -0.1681542 自分で相互相関係数を計算するときに必要なもの。 > mean.x <- mean(x) # x の平均（全部使う） > mean.y <- mean(y) # y の平均（全部使う） > ss.x <- var(x)*(n-1) # x の変動（全部使う） > ss.y <- var(y)*(n-1) # y の変動（全部使う） 相互相関係数を計算する。lag = 2~ 使用するデータベクトル~ x[(lag+1):n]-mean.x つまり，lag+1 から n まで，長さ n-lag。平均値 mean.x を引いておく。~ y[1:(n-lag)]-mean.y つまり，1 から n-lag まで，長さ n-lag。平均値 mean.y を引いておく。~ 両者掛け合わせて合計し（部分の共変動だ） sum((x[(lag+1):n]-mean.x) * (y[1:(n-lag)]-mean.y) ~ x，y 両方の全体の変動の積の平方根 sqrt(ss.x * ss.y) で割る > sum((x[(lag+1):n]-mean.x) * (y[1:(n-lag)]-mean.y) / sqrt(ss.x * ss.y) [1] -0.1681542 上の ccf の lag = 2 の下に書かれている数値 -0.1681542 と一致しているでしょ。~ つまり，ヘルプに書いている 「x[t+k] and y[t]」というのはそいういうこと，です。よ。~ x[(lag+1):n] と y[1:(n-lag)] どちらを基準にしていますか？~ そんなこと，どっちだっていいじゃないか。 -- &new{2020-03-13 (金) 11:11:38}; #comment **R のアップデート [#rac7c1ae] >[[私は誰でしょう？]] (2020-03-01 (日) 21:20:52)~ ~ 今まで，最新の R にするために（新しもの好きなので）~ あ，見ての通り，Mac ユーザです。~ ~ curl https://mac.r-project.org/el-capitan/R-3.6-branch/R-3.6-branch-el-capitan-sa-x86_64.tar.gz | sudo tar fvxz - -C /~ ~ ~ をやっていたのですが，2 ヶ月ほど前から，実行結果が思わしくなくなっています。~ つまるところ，件の URL はアクセス不可ということのようで，アップデートは出来ません。~ ~ 念の為~ http://mac.r-project.org/el-capitan/R-3.6-branch/R-3.6-branch-el-capitan-signed.pkg~ を使えばアップデートは出来ます。~ // - サイトの方のWAFとかでブロックされているのか, 単なる設定ミスかはわかりませんが当方でも同じ症状です. curl -O http://mac.r-project.org/el-capitan/R-3.6-branch/R-3.6-branch-el-capitan-signed.pkg sudo installer -pkg ./R-3.6-branch-el-capitan-signed.pkg -target / とかでお茶を濁しては... -- [[何処の誰かは知らないけれど]] &new{2020-03-02 (月) 16:03:50}; - 有難うございます。取りあえずは，それで，しのいでみます。 -- &new{2020-03-02 (月) 22:12:06}; #comment **連番でファイルを読み込みたいときにファイルに欠番があるときの処理 [#me7fd7f2] >[[Rビギナー]] (2020-02-26 (水) 10:12:45)~ ~ お世話になります．~ 下記のように1～3まで連番のファイルを読み込みたいのですが，例えば2番のファイルが無い場合にそのまま1，3のファイルを読み込んでプログラムが止まらないようにするにはどうしたらいいでしょうか．~ try関数で可能かと思いやってみたのですが，だめでした．~ ご教示ください．~ よろしくお願いいたします．~ ~ ＝＝ここから~ point <- 1:3 # for文を回す数を設定 for (i in point) { mypath <- file.path(getwd(), paste0(i,"t2-1.csv")) df <- assign(paste0("data",i), read.csv(mypath, header = T)) df <- data.frame(df) df head(df) df$X <- strptime(df$X,"%Y/%m/%d %H:%M'%S") head(df) df$X <- as.POSIXlt(df$X,"%Y/%m/%d %H:%M:%S",tz="GMT0") head(df) #https://qiita.com/tktz/items/733c37b1d6102ae52120 library(xts) df$X <- align.time(df$X-10*60/2,10*60) head(df) date = paste0("date",i) Temp = paste0("Temp",i) RH = paste0("RH",i) colnames(df) <- c(date,Temp,RH) head(df) mypath <- file.path(getwd(), paste0(i,"t1.csv")) write.csv(df,mypath) } // - insert 1 sentence. -- &new{2020-02-26 (水) 10:46:10}; mypath <- file.path(getwd(), paste0(i,"t2-1.csv")) if (! file.exists(mypath)) next - or -- &new{2020-02-26 (水) 11:09:20}; files = list.files(getwd(), pattern="[0-9]*t2-1.csv") for (f in files) { # describe what you want to do } - お二人の方，ありがとうございました．おかげさまでうまくいきました！ -- [[Rビギナー]] &new{2020-02-26 (水) 14:57:32}; - お二人の方，ありがとうございました．おかげさまでうまくいきました！ -- [[Rビギナー]] &new{2020-02-26 (水) 15:47:15}; #comment **時系列 [#h477736a] >[[マリン]] (2020-02-23 (日) 10:13:04)~ ~ ２つの時系列のデータ（約８０個）の相関（どのくらい類似性があるか）を~ 見たいのですが、方法を教えてください。~ // - まずはググってみて，それからにしなさい -- &new{2020-02-23 (日) 21:09:21}; - ？ -- [[マリン]] &new{2020-02-24 (月) 07:45:37}; #comment **固定効果を用いた回帰線 [#s967c99e] > (2020-02-14 (金) 14:34:30)~ ~ Book1というcsvファイルのデータを用いてGLMMによる解析を行ったのですが、固定効果の値を用いた回帰線がうまく引けません。~ csvファイルおよび解析を行ったスクリプトをBook1GLMMというファイル名でアップロードしましたので、ご教授いただけますと幸いです。~ よろしくお願いいたします。~ // - do! just shown as article '2020-01-30 (木) 12:37:42' the same program isn't it? -- &new{2020-02-14 (金) 17:11:14}; data <- read.csv("Book1.csv") data # update.packages # install.packages("lme4") library(lme4) par(family = "HiraKakuProN-W3") model1 <- glmer(POMN ~ log(siltclayN) + (1 | Place), family = Gamma(link = "log"), data = data) summary(model1) place = c("HA", "TN", "OG", "YN", "KG", "WA", "AS") # データフレームではこちらを使っているのでこれでなきゃだめ place2 = c("半田山", "田野", "小川", "与那", "上賀茂", "西粟倉", "足寄") # 凡例用には別のものを用意する pch = c(25, 18, 17, 16, 15, 9, 8) col = c("black", "red", "cadetblue", "magenta", "brown", "coral", "darkorchid") plot(data$siltclayN, data$POMN, pch = rep(pch, each = 5), col = rep(col, each = 5), cex = 3, cex.axis = 2) legend("bottomright", legend = place2, pch = pch, col = col, cex = 2) xx <- seq(min(data$siltclayN), max(data$siltclayN), length = 100) # yy <- exp(1.1136 + 0.4651 * xx) # そもそも間違い # lines(xx, yy, lwd = 2) for (i in seq_along(place)) { # 7 本の曲線を引くんでしょ？ newdata = data.frame(siltclayN=xx, Place=rep(place[i], 100)) yy = predict(model1, newdata=newdata, type="response") lines(xx, yy, col=col[i]) } - 7本の曲線は引けました。勉強不足でお恥ずかしいのですが、固定効果の推定値を用いて一本だけ回帰線を引くことはできないのでしょうか？ -- &new{2020-02-14 (金) 18:08:12}; - quit -- &new{2020-02-14 (金) 21:04:44}; #comment **ガンマ分布 [#f464b927] > (2020-02-06 (木) 21:14:35)~ ~ GLMMなどで、目的変数に用いるデータが連続値で0以上であればガンマ分布を用いると統計の教科書等には書いてあるのですが、連続値で0以上(データの最小値は0)のデータでガンマ分布を用いて解析しようとすると~ 「non-positive values not allowed for the 'gamma' family」というエラーが出ます。~ また、同じデータで線形混合効果モデルなら結果は出るのですが、「boundary (singular) fit: see ?isSingular」という警告？が出るせいなのかランダム効果を考慮した回帰線が一本しか出ません。~ 解決策などありましたら、お教えいただけますと幸いです。~ // - show your data and your script(program). -- &new{2020-02-06 (木) 22:07:07}; - model1<-glmer(A~log(B)+(1|Place),family=Gamma(link="log"),data=d) In this script was error 「non-positive values not allowed for the 'gamma' family」. -- &new{2020-02-06 (木) 22:26:52}; - model2<-lmer(A~B+(1|Place),data=d) In this script was 「boundary (singular) fit: see ?isSingular」 -- &new{2020-02-06 (木) 22:28:59}; - The data was uploaded. It is data named data.csv. -- &new{2020-02-06 (木) 22:31:00}; - 'non-positive values not allowed for the 'gamma' family' means 0 is not non-positive, isn't it. > library(lme4) > d = read.csv("data.csv") > model1 <- glmer(A ~ log(B) + (1|Place),f amily=Gamma(link="log"),data=d) eval(family$initialize, rho) でエラー: non-positive values not allowed for the 'gamma' family using d2 > d2 = d[d$A > 0, ] > model1-2 <- glmer(A~log(B) + (1|Place), family=Gamma(link="log"), data=d2) > model1-2 no error happend.&br; can you understand?&br;you must *respect* 'error message'. it warns you are wrong, and he says the reason 'why'.&br;if your data contains 'zero's, it may not be a gamma distribution!!!!! -- &new{2020-02-06 (木) 23:37:25}; - I will try it once. One more thing I want to ask. -- &new{2020-02-07 (金) 00:00:48}; - Now I uploaded another data. When the script of model2 is executed with this data, the message "boundary (singular) fit: see? IsSingular" appears, and only one regression line considering random effects is drawn. What is the cause? -- &new{2020-02-07 (金) 00:03:40}; - you! just , imput '? IsSingular'(without last '?') in console!!!. message says so, ah. -- &new{2020-02-07 (金) 00:11:48}; model2<-lmer(NitMAOMrecover~MAOMCN+(1|Place),data=example) summary(model2) place=c("HA","TN","OG","YN","KG","WA","AS") pch=c(25,18,17,16,15,9,8) col=c("black","red","cadetblue","magenta","brown","coral","darkorchid") plot(example$B,example$A,pch=rep(pch,each=5),col=rep(col,each=5),cex=3,cex.axis=2) legend("bottomright",legend=place,pch=pch,col=col,cex=1.5) xx<-seq(min(example$B),max(example$B),length=1000) for(i in seq_along(place)){ yy<-predict(model2,newdata=data.frame(B=xx,Place=rep(place[i],1000)),type="response") lines(xx,yy,col=col[i],lwd=2)} - Excuse me. I want to show the script to draw the regression line on the site, but I don't know how to do it. -- &new{2020-02-07 (金) 00:29:15}; - they said. "その前に，「[[投稿における注意事項]]」を読んでから！" and "投稿する前にまず &heart; [[投稿文書の書式:http://www.okadajp.org/RWiki/?%E6%95%B4%E5%BD%A2%E3%83%AB%E3%83%BC%E3%83%AB]] &heart; を読んでください" -- &new{2020-02-07 (金) 08:20:14}; - surprisingly, it's a reckless and messy program.~ Four variables "Sample,Place,A,B" are included in your "example.csv".~ but, your wrote "NitMAOMrecover~MAOMCN+(1|Place),data=example".~ more over, you wrote variables "A", "B" as "example$B" or "example$A"~ a data-frame in lmer function is different from a "new" data-frame in predict function. ~ you have to study a basis of R. -- &new{2020-02-07 (金) 08:23:45}; #comment **GLMMについて [#tc8ba31b] >[[ra]] (2020-01-27 (月) 11:58:45)~ ~ R初心者です。~ GLMMで解析を行った後、散布図に回帰線を引きたいのですが、predict関数を読み込ませようとすると（オブジェクト"Place"がありません）というエラーが出ます（Placeはランダム効果です）。このエラーを解消し、回帰線を引くにはどうすれば良いでしょうか？~ 入力したコマンドを載せておきます。~ よろしくお願いいたします。~ ~ model1<-glmer(POMN~log(siltclayN)+(1|Place),family=Gamma(link="log"),data=data)~ plot(data$siltclayN,data$POMN,pch=c(25,18,17,16,15,9,8),cex=3,cex.axis=2)~ legend("bottomright",legend=c("HA","TN","OG","YN","KG","WA","AS"),pch=c(25,18,17,16,15,9,8),cex=2)~ xx<-seq(min(data$siltclayN),max(data$siltclayN),length=100)~ yy<-predict(model1,newdata=data.frame(siltclayN=xx),type="response")~ lines(xx,yy,lwd=2)~ // - new data set has to include 'siltclayN' *and* 'Place' -- &new{2020-01-27 (月) 13:47:56}; - Sorry,how do I create that data set? -- &new{2020-01-27 (月) 16:26:31}; - just, create -- &new{2020-01-27 (月) 22:11:04}; - POMN <- c(data$POMN) siltclayN <- c(data$siltclayN + Place) dataframe <- data.frame(POMN, siltclayN) model1 <- glmer(POMN~log(siltclayN) + (1|Place1), family=Gamma(link="log"), data=dataframe) plot(dataframe$siltclayN, dataframe$POMN, pch=c(25,18,17,16,15,9,8), cex=3, cex.axis=2) legend("bottomright", legend=c("HA","TN","OG","YN","KG","WA","AS"), pch=c(25,18,17,16,15,9,8),cex=2) xx <- seq(min(dataframe$siltclayN), max(dataframe$siltclayN), length=100) yy <- predict(model1, newdata=data.frame(siltclayN=xx), type="response") -- &new{2020-01-27 (月) 23:14:23}; - What's wrong? -- &new{2020-01-27 (月) 23:16:01}; - newdata=data.frame(siltclayN=xx, Place = ???????) -- &new{2020-01-28 (火) 10:50:09}; - In this command and Place=xx, an error has occurred. Where should I improve? -- &new{2020-01-28 (火) 14:42:23}; - Place is factor isn't it? 'Place = xx' may be invalid!&br; I'll show you an e xample. > df = data.frame(x=c(3,2,1,3,5), f=factor(c(1,2,2,3,3)), y=c(2,4,7,5,6)) > fit = lm(y ~ x + f, data=df) # f may not be a factor > predict(fit) # this prediction uses data-frame df 1 2 3 4 5 2.0 5.4 5.6 5.7 5.3 if you predict with another data-frame, the data-frame must contain 'x' and 'f' following newdata contains 'x' and 'f'. first 5 elements are same in df(df$x, df$f), last 3 elements is exactly new data > newdata = data.frame(x=c(3,2,1,3,5, 2,3,5), f=factor(c(1,2,2,3,3, 1,2,2))) > predict(fit, newdata) 1 2 3 4 5 6 7 8 2.0 5.4 5.6 5.7 5.3 2.2 5.2 4.8 # first five values are identical to previous prediction. last 3 values are new prediction. do you understand?? -- &new{2020-01-28 (火) 16:28:06}; - Place is a factor. In that case, is Place=c(data$Place) ok? Or should all numerical data and factor data be written in c of x,f and y? -- &new{2020-01-29 (水) 10:40:07}; - I ran it with reference to your example, how can I draw a prediction line from here? -- &new{2020-01-29 (水) 11:37:47}; - if new data-frame is declared as 'newdata=data.frame(siltclayN=xx, Place=pp), pp[i] corresponds to xx[i].&br;examine above example carefullly. predict function must be called as 'predict(fit, newdata=data.frame(x=c(2,3,5), f=c(1,2,2)))' -- &new{2020-01-29 (水) 13:34:36}; - 'prediction line'?? no, no, prediction *prediction plane*. there are *two* !! independent variables(x, f) and one dependent variable(y), so, it's a 3 dimensional data(Actually, if f is a factor and has more than 2 levels, the dimension is more than 3!!). your data is 3 or more dimensinal, too.&br;you know, in case of simple regression analysis y = a*x+b you have one independent variable x, and one dependent variable y, it's a *two* dimensional data, and you can draw two dimensional scatter plot or regression *line*. &br;if you have 3 independent variable model, no one can draw predicted hyperplane(you can't imagine 4 dimensional space, can you??.&br;however, if you draw POMN vs. siltclayN by Place, you can draw 2 dimensional scatter or regression line. -- &new{2020-01-29 (水) 13:41:39}; - your model is POMN = α * log(siltclayN[i]) + β * Place[i]. if you want to predict POMN for new data set, you must prepare new data pair siltclayN[j] and Place[j], j = 1, 2, 3, ... -- &new{2020-01-29 (水) 14:02:09}; - your data-frame is curious, isn't it?&br; what do you mean siltclayN <- c(data$siltclayN + Place) dataframe <- data.frame(POMN, siltclayN) if Place is a factor variable you cant add to data$siltclayN. 'data$siltclayN + Place' cause an error. -- &new{2020-01-29 (水) 14:47:36}; - upload your data is best way to get solution -- &new{2020-01-29 (水) 14:54:46}; - Did you fall into a deadlock? -- &new{2020-01-29 (水) 22:32:29}; - yes... I don't know how to draw a regression line considering random effects. -- &new{2020-01-30 (木) 10:34:08}; - Can you tell me what to improve specifically? -- &new{2020-01-30 (木) 10:36:21}; - upload your data is best way to get solution -- &new{2020-01-30 (木) 10:46:16}; - Upload now. It is a file called Book1. -- &new{2020-01-30 (木) 10:49:42}; - alright, I was able to upload. -- &new{2020-01-30 (木) 11:31:58}; - Your program is strange overall. for example, pch=c(25,18,17,16,15,9,8) must be pch=rep(c(25,18,17,16,15,9,8), 5). examine each point's coordinates.&br; you may want a result like this. -- &new{2020-01-30 (木) 12:37:42}; library(lme4) data = read.csv("Book1.csv") model1 <- glmer(POMN ~ log(siltclayN) + (1 | Place), family = Gamma(link = "log"), data = data) summary(model1) place = c("HA", "TN", "OG", "YN", "KG", "WA", "AS") pch = c(25, 18, 17, 16, 15, 9, 8) col = c("black", "red", "cadetblue", "magenta", "brown", "coral", "darkorchid") plot(data$siltclayN, data$POMN, pch = rep(pch, each = 5), col = rep(col, each = 5)) # , cex=3, cex.axis=2) legend("bottomright", legend = place, pch = pch, col = col) # ,cex=2) xx <- seq(min(data$siltclayN), max(data$siltclayN), length = 100) for (i in seq_along(place)) { yy <- predict(model1, newdata = data.frame(siltclayN = xx, Place = rep(place[i], 100)), type = "response") lines(xx, yy, col = col[i]) } &ref(glmer.png); - I wanted to make this diagram. It was very helpful! I didn't have enough knowledge. -- &new{2020-01-30 (木) 13:48:53}; - I'm really thankful to you! -- &new{2020-01-30 (木) 13:50:02}; - I want to draw a regression line with the estimated fixed effects of model1, but I cannot do it well. -- &new{2020-02-13 (木) 21:49:52}; #comment **文字列からバイト数を指定して切り出す [#dcc3e621] > (2020-01-26 (日) 22:45:01)~ ~ 半角・全角が混在する文字列からバイト数を指定して、部分文字列を抽出したいのですが、どのような関数がありますでしょうか。~ 例えば~ x <- "あaいbうc"~ という文字列から3-5byte目を抽出したいです。~ （"aい"という部分文字列を抽出したいです。）~ ~ substringは文字数を指定する関数のため、byte単位で抽出が不可能なようです。~ ~ お手数ですがご教示お願いいたします。~ // - ? substr / If an input element has declared "bytes" encoding (see Encoding, the subsetting is done in units of bytes not characters. -- &new{2020-01-27 (月) 11:46:10}; - そもそも3-5byte目が「aい」というは正しいですか？ -- [[S]] &new{2020-01-27 (月) 16:29:39}; > x <- "あaいbうc" > charToRaw(x) [1] e3 81 82 61 e3 81 84 62 e3 81 86 63 「あ」は3バイト文字なので途中で切られて、3-5byte目は「82 61 e3」です。これは「あ」の3バイト目、1バイト文字の「a」、「い」の1バイト目です。~ 全角半角の字数をバイト数で何とかしようとするのは無理があ流のではないでしょうか。 - it's just he wants! -- &new{2020-01-27 (月) 22:13:13}; - euc-jp uses 2 bytes for 1 japanese character > x <- "あaいbうc" # UTF-8: あ == \xe3\x81\x82 > xx = iconv(x, "utf-8", "euc-jp") > xx [1] "\\xa4\\xa2a\\xa4\\xa4b\\xa4\\xa6c" # EUC-JP: あ == \xa4\xa2 > Encoding(xx) = "bytes" > iconv(substr(xx, 3, 5), "euc-jp", "utf-8") [1] "aい" cp932 uses 2 bytes for 1 japanese character, too > x <- "あaいbうc" > yy = iconv(x, "utf-8", "cp932") > yy [1] "\x82\xa0a\x82\xa2b\x82\xa4c" # CP932: あ == \x82\xa0 > Encoding(yy) = "bytes" > iconv(substr(yy, 3, 5), "cp932", "utf-8") [1] "aい" do you understand? -- &new{2020-01-28 (火) 11:38:51}; - 文字コードを変換してsubstrで切り出すということですね。ご回答ありがとうござました。 -- &new{2020-01-30 (木) 00:44:16}; #comment **apply系関数でfunに線形回帰lmが適用可能な関数 [#v28d972a] >[[TY]] (2020-01-07 (火) 15:59:04)~ ~ シミュレーションで同じ関数を繰り返し(10^5~10^6回程度)実行する場面を想定しています．~ 例えば，以下のsample_dataでREP毎にlm(y~x)を実行したいです．~ for文で実行することはできるのですが，一般にapply系関数の方が実行が速いと聞き，~ 探しています．apply系関数に関する知識は，tapply関数にmean関数を適用できる，~ tapply(sample_data$y, sample_data$REP,mean)~ 程度のレベルです．~ ~ 【sample_data】~ y x REP 1 133 3.9 1~ 2 81 -2.3 1~ 3 104 -1.2 1~ 4 89 0.7 1~ 5 142 4.0 1~ 6 74 -2.9 1~ 7 140 3.9 1~ 8 170 4.4 1~ 9 128 1.6 1~ 10 107 1.2 1~ 11 49 -4.3 2~ 12 79 -2.9 2~ 13 77 -3.2 2~ 14 125 1.8 2~ 15 82 -1.1 2~ 16 132 2.6 2~ 17 99 -0.0 2~ 18 140 2.1 2~ 19 140 4.9 2~ 20 88 -1.1 2~ // - 「REP毎にlm(y~x)を実行したい」のか，同じような 10^5 ほどのデータを分析したいのか。教えてあげてもいいけど，そのへんどうなのよ。~ つまり，rep が 1~10^5 くらいの種類があって，それぞれについて x,y の組のデータ（サンプルサイズが 10 くらい）があるんだけど，それぞれについて lm で分析してータセットごとの結果が欲しいのか？&br;結果を待っているより，やってみた方が速いよ（ワラ）。反応があるという保証がない限り，回答しない。&br;それにしても，質問の仕方がへた。 -- &new{2020-01-10 (金) 22:10:35}; - ありがとうございます．出直してきます． -- [[TY]] &new{2020-01-13 (月) 10:59:17}; - 必要なのは，そういう対応ではないと思うが。質問にも答えないし，誠意がない。というか，本当に回答が必要なのか？ -- &new{2020-01-16 (木) 21:03:50}; #comment
