COLOR(blue){SIZE(20){R および RjpWiki に関する質問コーナー}}~
COLOR(red){SIZE(20){このコーナーはスパム対策のため閉鎖します。以後は [[Ｑ＆Ａ(中級者コース):http://www.okada.jp.org/RWiki/index.php?%A3%D1%A1%F5%A3%C1%20%28%C3%E6%B5%E9%BC%D4%A5%B3%A1%BC%A5%B9%29]]をご利用ください}}~
//COLOR(red){注意}：[[Ｑ＆Ａ]] コーナーが少し重くなり過ぎた気がしますので、Ｑ＆Ａ(2) として続けたいと思います（画像ファイルと一緒に引越しするのが困難なので、通し番号方式でいきたいと思います）。~
COLOR(red){注意}：新規記事用の入力欄は以下の目次の直後にあります。各記事への追加コメントの入力欄は、各記事の後にあります。質問が長くなる場合は、まず出だしだけ投稿し、次に編集ボタンでそれを修正、追加するのが便利でしょう。
もし新規に質問をされる場合は、効率的にコメントを得るためにも、善意のコメント者の時間を無駄にしないためにも、以下に特にご注意下さい。[[R のメイリングリストへの投稿記事のガイド]] も参考になるでしょう。
- 広い意味で R に関係する話題に限定!
- 学校の宿題は自分で考えるべきです
- 具体的な状況がわかる様に背景説明をけちらない
- [[CRAN の検索エンジン:http://finzi.psych.upenn.edu/search.html]] や [[マニュアル:http://buran.u-gakugei.ac.jp/~mori/LEARN/R/]] を調べた上での質問かどうかを書き添えれば、空振りコメントが減るでしょう。また必須ではありませんが r-help の投稿ガイダンス [[R のメイリングリストへの投稿記事のガイド]] も参考にして下さい。
- 自分には馴染みでも、他の人には馴染みの無い手法・概念もあります。日本語は得てしてわかりにくい(英語名をそえると案外わかるものです)ものです。また英語省略名は便利ですがしばしば仲間うちだけでしか通用しないと考え、一般的な言い替え、説明も考えて下さい。
- もし問題が解決した(結局解決しない場合を含め)場合はその旨報告していただくとコメント者や他の人に参考になります。また後で関連情報をこのWikiに要約記事としてご投稿いただくと今後の参考になります。
- 後で他の人が参照する際の便宜のため、タイトルには簡潔で内容が推察できるようなものをつけて下さい~
- 回答者は神通力を持ってはいません。「これこれしたが、うまく動きません」といった質問には、あなたの試みたことが回答者に再現できるような、データ、コード、エラー出力(そして、R のバージョン、使用 OS、使用非標準パッケージ名)を添えることをお勧めします。あまりに長くなる場合は適宜編集して下さい。~
~
----
アーカイブ：[[Ｑ＆Ａ(旧1)]]~
COLOR(red){SIZE(20){このコーナーはスパム対策のため閉鎖します。以後は [[Ｑ＆Ａ(中級者コース):http://www.okada.jp.org/RWiki/index.php?%A3%D1%A1%F5%A3%C1%20%28%C3%E6%B5%E9%BC%D4%A5%B3%A1%BC%A5%B9%29]]をご利用ください}}~
#contents
----
COLOR(red){注意}：投稿失敗した記事をページ先頭の「新規｜編集｜差分｜添付」の「編集」から消すときは、#article と書いた行まで消さないように注意しましょう。新規の記事投稿ができなくなります。
#article
//#article
**文字式による数式の関数認識
>[[KOSHI]] (2006-01-18 (水) 22:58:51)~
~
expression()の形式で入力した数式を，関数として認識したい場合，~
例えば
> f <- expression(x*y)
を（x,y）の関数に直したいとき，
> F <- function(x,y){}
> body(F) <- f
とやれば
> F <- (x, y) x * y
とできます．しかし，これはfで使われる変数が（x,y）と既知であるから出来る操作です．~
もし，fでの変数が未知の場合，これを関数認識するにはどうしたらよいでしょうか？つまり，ある関数Gで，
> g <- G(expression(x*y*z*u))
などとすると,ちゃんとこの4つの変数（x,y,z,u）を使って~
> g <- function(x,y,z,u) x*y*z*u
となるような関数を作ることは可能でしょうか．~
~
ちなみに~
deriv(expr, name, func=T)
という関数を使うと，勝手に書いた数式~
expr=expression(...)
の，変数nameによる1階偏微分を，その文字による関数として認識させているようです．こういうことを自分でもやりたいのですが，derivのソースを見ることができず，結局わかりません．~
~
どなたかお力をお貸しいただけないでしょうか．~
//
-ソースを見ることについて，どこかのページにも書いてありますが，底に書いてあることもちょっと違うこともあるかもしれないので，ここにまとめてみようかなぁ~
まず，関数名だけをタイプしてみる
> deriv
function (expr, ...)
UseMethod("deriv")
<environment: namespace:stats>
なんだぁ？これだけか？と思ったら，
> methods(deriv)
[1] deriv.default deriv.formula
と，やってみる。ここに示されるのが，実際に呼び出される候補の関数名。~
単に関数名をコンソールにタイプするのではなくて，以前に出力されて namespase: の次にある名前をまず書いて，その後にコロン（：）を三つ書いてその後に，その関数名を書く
> stats:::deriv.default
function (expr, namevec, function.arg = NULL, tag = ".expr",
hessian = FALSE, ...)
.Internal(deriv.default(expr, namevec, function.arg, tag, hessian))
<environment: namespace:stats>
もう一つ
> stats:::deriv.formula
function (expr, namevec, function.arg = NULL, tag = ".expr",
hessian = FALSE, ...)
{
if ((le <- length(expr)) > 1)
.Internal(deriv.default(expr[[le]], namevec, function.arg,
tag, hessian))
else stop("invalid formula in deriv")
}
<environment: namespace:stats>
これらの関数のどれがどのように使われるかは，help で検索のこと。~
ここで示されるのは，R の範囲で書かれるプログラムであって，複雑な計算処理プログラムは更にここから呼ばれる。上の例では .Internal(deriv.defarult.... などとある部分かな。~
そのような部分についてのプログラムの詳細は，（えい面倒だなぁ。。），普通のダウンロード・インストールでは実装できないことに依存するので，，，というか，そこまで望む人には，ソースを探索する必要がありますということだけを言えばいいだろうし，そういわれても分からない人は，ソースをみても特に情報は得られないのではないか（ご，ごめんなせい，，，石を投げないでおくんなせーーー）。~
てなことかな（はずしたか。。） -- &new{2006-01-18 (水) 23:16:34};
-それとも，あなたは，このようなことがお望みなのかな？ as.function
> func <- as.function(a <- function(x, y) x+y)
> func
function(x, y) x+y
> func(2,4)
[1] 6
> func <- as.function(a <- function(x, y) x*y)
> func(3,7)
[1] 21
一応そういうことで。 -- &new{2006-01-18 (水) 23:36:17};
-私の場合，ソースを見ても分からない部類だと思います・・・．as.functionも，結局はどの変数を使うかを知っていて定義しているのだと思います． 私が作りたいのは，ユーザーが勝手に作った数式を，関数に変換するプログラムで，たとえユーザーが，どんな文字を使おうと，それを一発で関数にしてしまうようなものです． -- [[KOSHI]] &new{2006-01-18 (水) 23:40:50};
-実際，deriv（）はどんな文字を用いて数式を書いても，それをその文字で微分して，その文字の関数に出来ています．なんででしょう・・・？ -- [[KOSHI]] &new{2006-01-18 (水) 23:47:02};
-具体的には，どんなものですか？上の，func だって，x,y を使って定義されようが，a,b を使って定義されようが同じように動きますが。う〜〜っん。ユーザが"a*x^2+b*x+c", という文字列を入れて，x が変数，a,b,c が定数，そして，関数名は quad にしたいなあといったら，quad <- function(x) a*x^2+b*x+c という関数定義がしたいなぁということですか？そのようなことであるとしても，特定の限定された条件下ではいろいろ実現方法はあるかもしれませんしねえ。 -- &new{2006-01-18 (水) 23:48:42};
-＞実際，deriv（）はどんな文字を用いて数式を書いても，それをその文字で微分して，その文字の関数に出来ています．なんででしょう~
derive を読んでください(^_^;)~
もしくは，ウイザード召喚！！(^_^;) -- &new{2006-01-18 (水) 23:49:53};
-たとえば，ユーザが f <- expression(x*y) と入力することが分かっていたら，
> G <- function(expr){
F <- function(x,y){ }
body(F) <- f
return(F)
}
というGを用意しておけば(x,yの数式なら何でも)関数に出来ますが，このGを用意したところで，ユーザが
> f <- expression(u*v)
と書いてしまえば，G(f)は(u,v)の関数に出来ないじゃないですか．ここでG(f)がちゃんと(u,v)の関数にもなって欲しいんです．そういうGをつくりたいんですが・・・こんな感じで意味分かりますでしょうか・・・？ -- [[KOSHI]] &new{2006-01-19 (木) 00:01:06};
-関数定義に出てくる名前をパースして，それぞれが想定しているどれにあたるか判定して，，，言うはやすく，実現するは難しですか。。。というか，ユーザがx*yでなくてx+yと書いたら，更に x+y+z と書いたらとか，そういうことはなしですか。。。そういうこともありなら，やはり，もとの C プログラム（かどうか，知りませんし，知りたいとも思いませんが）その方面を攻めないといけないかもしれませんね。 -- &new{2006-01-19 (木) 00:08:20};
-deriv も何が変数（とみなすべき）かを指示する仕様になっていますよね。汎用性のある手順を考えているならこれは避けられないはず。ということは、以下のようなやりかたが望みうる限界なのでは？ alist を文字列から合成する方法はわからぬ。-- &new{2006-01-19 (木) 16:49:08};
> a <- 3; b <- 5
> f <- alist(x=,y=a,x*y+b)
> F <- as.function(f)
> F
function (x, y = a)
x * y + b
> F(1)
[1] 8
-＞ alist を文字列から合成する方法.....一時ファイルに書き出して，source で読む。。バキッ;_; -- &new{2006-01-19 (木) 17:36:10};
-なるほど。一時ファイル名をその場で作る関数が確かあったはず。 --?{2006-01-19 (木) 20:56:34};
-tempfile 関数ですね。。 -- &new{2006-01-20 (金) 10:56:42};
-皆様、いろいろとご助言有難うございます。自分でも少し進展して、ユーザが入力する数式の変数の個数を予め知っていれば、derivを利用して、たとえユーザがどんな変数でどんな数式を入力しようと、それを関数化することに成功しました。例えば、ユーザが２個の変数を持つ数式しか入力しないと分っている時、
> Func <- function(expr){
df <- deriv(expr, all.vars(expr), func=T) #とりあえず全ての変数で微分しておく。
F <- function(p,q){ #ここを変数２個にしておく。
df(p,q)[1] #derivの[1]には微分しない関数が返されていることを利用。
}
return(F)
}
こうして関数Fを作れば、ユーザが入力する「変数の個数」さえ分っていれば、こっちの好きな変数で関数化することが出来るように思います。~
あとは、自然数nを与えられた時に、自動的にn変数の空の関数を作るような関数（すなわち、G(n)=function(x1,...xn）{}となるような関数）が構成出来れば、n <- length(all.var(expr))とすることによって、
> body(G(n)) <- df(x1,...xn)[1]
で完成するのではないでしょうか？~
しかし、こういうG(n)ってできますでしょうか？~
forとか使っても関数の変数を増やせなくて・・・あと一歩だと思うのですが。~
よろしくお願いします。 -- [[KOSHI]] &new{2006-01-20 (金) 18:59:21};
-上に示唆されていた方法を使えば（ファイルはテキスト、文字列、だという特徴を使う）。それと編集の際は最後の #comment を残しておきましょう。 -- [[間瀬茂]] &new{2006-01-21 (土) 11:24:12};~
> G <- function(n) {
tmpfile <- tempfile() # 一時ファイル名を作る
.x <- paste("x", 1:n, ",", sep="") # 引数リスト文字列ベクトル
.x[n] <- substr(.x[n], 1, nchar(.x[n])-1) # 最後の引数文字列からカンマを取る
sink(tmpfile)
cat("function (", .x, ") {}?n") # 関数定義文字列をファイルに出力
sink() # ファイル出力を閉じる
.f <- source(tmpfile)[[1]] # ファイルから関数定義を読み込む
unlink(tmpfile) # 一時ファイルを削除
return(.f)
}
> f5 <- G(5)
> f5
function ( x1, x2, x3, x4, x5 ) {}
> str(f5)
function (x1, x2, x3, x4, x5)
- attr(*, "source")= chr "function ( x1, x2, x3, x4, x5 ) {}"
- G<-function(n)eval(parse(text=paste("function(",paste("x",1:n,sep="",collapse=","),"){}"))) -- [[takahashi]] &new{2006-01-21 (土) 18:53:25};
-↑ 上手い！ しかし、ファイルに一旦書き出す方法は（少し野暮ったくみえるけれど）汎用性がある、覚えておいて良い方法ですね。 -- [[間瀬茂]] &new{2006-01-21 (土) 23:18:18};
-どうも，このページに粘着質のスパム投稿者が住み着いたようで。#comment も意図的に全部はずす羽目になって。そうすると今度はスレッド立てるし。どうしょうもないね。 -- &new{2006-01-23 (月) 09:01:41};
-Q & A (中級) とでも改名して再出発しましょうか。みたところ過去記事の AIDS というキーワードでどこかの宛名業者に覚えられてしまったようですね。 -- &new{2006-01-23 (月) 09:43:06};
-スパム対策のためＱ＆Ａコーナーは閉鎖します。以後は [[Ｑ＆Ａ(中級者コース):http://www.okada.jp.org/RWiki/index.php?%A3%D1%A1%F5%A3%C1%20%28%C3%E6%B5%E9%BC%D4%A5%B3%A1%BC%A5%B9%29]] をご利用ください。 -- [[間瀬茂]] &new{2006-01-23 (月) 12:57:52};
**共分散分析
>[[S.O]] (2006-01-13 (金) 03:10:39)~
~
共分散分析のやり方が、いまいち分かりません。~
yが目的変数で、x1,x2は連続量で、x3,x4は離散です。~
共分散分析は、lmよりaovの方がよいと言われていたので、以下の式で共分散分析ができているのでしょうか？~
aov(y~x1+x2+x3+x4)~
この場合、普通の４元配置分散分析だと思い、x1,x2を共変量にしないといけないのかと思ったのですが、やり方がよく分かりません。~
お手数ですが、どうぞ宜しくお願いいたします。~
//
-「いまいちわからない」ではこっちも「いまいちわからない」です。それとgoogleとかを使って調べたのでしょうか？いろいろ落ちているように思いましたが。 -- &new{2006-01-13 (金) 08:38:28};
**Mac版のR2WinBUGS
>[[まさひと]] (2006-01-11 (水) 20:54:49)~
~
パッケージではMac用のR2WinBUGSが配付されいますが、もともとWinBUGS自体がwindows用なのでうまく呼び出せないと思うのですがいかがでしょうか？ちなみに自分は呼び出しに失敗しました。~
//
-example(bugs) で，なんか動いているけど。よくわかりませんね。 -- &new{2006-01-11 (水) 23:20:23};
-wine経由で動かす(system関数使って)物です. 私の手元でも(Debianですが)wineでWindows版Rが動いてます. example(bugs)で動くのはパッケージの自動ビルドの為にドキュメント(Rd)に ?examples 内のコードに ?dontrun で囲むとコメントにされて実行されなくなります. wineは一時からみると随分進化しましたから, 是非お試しあれ. -- [[なかま]] &new{2006-01-12 (木) 16:13:19};
-Darwineをインストールしてみました. GUIのバージョンの関係で10.3は動作しませんでしたが, 10.4ではちゃんと動作しました. X11環境(元はwineですから)で動くので若干使いにくさは感じますが, それでも動きます(凄すぎる...). :-) -- [[なかま]] &new{2006-01-13 (金) 10:59:37};
-なかまさん、ありがとうございます。さっそく試してみます。 -- [[まさひと]] &new{2006-01-13 (金) 11:42:50};
-あー、すいません。よく見るとまだまだ普通には動かんようです。->darwine -- [[なかま]] &new{2006-01-13 (金) 15:14:29};
**Rで＠Risk的なパッケージを自作
>[[Tak]] (2006-01-10 (火) 18:43:58)~
~
リスク解析の勉強をかねて、Rで＠Riskの様なパッケージを自作しようと思っています。＠RiskはPalisade社が販売している、リスク解析の為のMS Excelアドインです。モンテカルロ・シュミレーション及び、sensitivity analysisやdecision tree analysisなどができるようです。~
http://www.palisade-europe.com/~
~
既存のRのパッケージをできるだけ使って、＠Riskの様な包括的なリスク解析パッケージが作れたらいいと思っています。手始めに、既存のRのパッケージやファンクションで、＠Risk的なリスク解析に役立ちそうな物があったらおしえていただきたいのです。例えばcor.testのmethod = 'spearman'はそのままnonparametric sensitivity analysisに使えるでしょうし、RExcelはExcelをインターフェイスにするのに役立ちそうです。~
~
他にparametric sensitivity analysis、モンテカルロ法、decision tree analysisなどなどのお勧めパッケージを教えてください。よろしくお願いします。~
//
-Palisade社が何らかの特許を申請していた場合、完全移植には注意が必要です。リスク解析であればPalisade社の名前をあげる必要はなかったのでは？ -- &new{2006-01-10 (火) 19:49:39};
-確かにPalisade社の名前をあげる必要はなかったですね。＠Riskの完全移植を狙っているわけではありません。ただ、モンテカルロ・シュミレーションやsensitivity analysisを誰でも簡単にできるフリーパッケージを作れたらいいなと思いました。 MS Wordに対抗したOpenOfficeの様に似てしまうかもしれませんが、今の段階では先走りになってしまいますね。こんなところでどうでしょう。-- [[Tak]] &new{2006-01-10 (火) 20:06:42};
**plotのwidthの種類
>[[浩太]] (2006-01-05 (木) 13:52:55)~
~
barplotの棒の幅を定めるwidthって何種類ぐらいの幅がありますか？~
具体的にどんな種類の幅がありますか？~
教えてください~
//
-質問の意味がわかりません。widthは数値ベクトル（リサイクルされる）で決められるものですから，種類は無限個あるとしか言えないでしょう。
> x <- c(1,4,8,5,2)
> barplot(x, width=c(1.2, 3.4, 5.321, 2.345))
これにより以下のような図が得られるんだが。。。
#ref(barplot.png)
それとも，
> barplot(x, width=1)
> barplot(x, width=1000)
のような描画では差がないように見えるからということかな？~
ヘルプを読めばわかるように，xlim で横軸の長さを決めない限り width による差は目に見えない。だって，横幅がウインドウに収まるように拡大・縮小されるわけであるから。~
Specifying a single value will no visible effect unless xlim is specified. -- &new{2006-01-05 (木) 14:35:56};
-ｐarを使って図を重ね合わせる時に差をつけるときも数値でやるんですか？ -- [[浩太]] &new{2006-01-05 (木) 15:23:36};
-浩太さんがparで図を重ねることができる人ならば、実際に試してみればよいのでは？その上の不具合であれば、「〜だが…できない」と具体的に質問した方が良いでしょう。既に幅の変更方法はわかったのですから。 -- &new{2006-01-05 (木) 15:40:39};
-四つのbarplotをparを使って重ねあわせました。それぞれの違いをはっきりさせるために色の違い以外に幅の違いをだすためにwidthを100,50,10,5にしたんですが違いがうまれませんでした。 -- [[浩太]] &new{2006-01-05 (木) 16:53:55};
-用いたプログラムをここに載せればいかがですか。長すぎたり，データが複雑なら，単純化してなおかつ状況を再現できるものを。~
なお，投稿の失敗は，上のほうの「編集」をクリックすることで編集できますよ。質問の仕方も練習して，上手な質問・後からこれを見る人の参考になるような質問にしましょう。 -- &new{2006-01-05 (木) 18:33:23};
-一般的に言えば，幅を変えると必然的に面積も変わるわけで，人間は同じ高さでも幅が広い方の長方形が表すものが大きいと受け取りますよ。あなたのやろうとしていることは「統計学で嘘をつく法」の実践になっているのではないかと危惧します。どうしてもやりたいのなら，止めはしませんが(^_^;) -- &new{2006-01-05 (木) 18:43:27};
**barplotでNameが表示されない
>[[鍋]] (2006-01-03 (火) 14:29:13)~
~
> mydata
[1] 29 24 20 18 25 22 23 20 20 23 20 30 21 20 19 26 14 23 17 16 13 10 9 10 7
[26] 2 2 4 1 3
> main
[1] "Turnover Number"
> names.arg
[1] 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
[26] 26 27 28 29 30 31 32
> xlab
[1] "Event Time Number"
> ylab
[1] "Turnover"
> ylim
[1] 0 32
> barplot(mydata, main=main, names.arg=names.arg, xlab=xlab, ylab=ylab, ylim=ylim)
Error in barplot.default(mydata, main = main, names.arg = names.arg, xlab = xlab, :
incorrect number of names
なぜこのようなErrorが表示されるのかわかりません。~
incorrectになるのはどういったときなのでしょうか？~
//
-名前の数が合わないと申しておりますし,実際合ってないように見えるんですがいかがでしょう? -- [[なかま]] &new{2006-01-03 (火) 15:34:49};
-length(mydata) の結果と，length(names.arg) の結果を比べてみましょう。名前の数が正しくないと言われているのですから，原因はわかるでしょうし，その対処法もわかると思いますが？？ -- &new{2006-01-03 (火) 15:35:56};
-国際化された R を使えば，エラーメッセージも日本語で出るので，よろしいのではないのでしょうかねぇ。
以下にエラーbarplot.default(mydata, main = main, names.arg = names.arg, xlab = xlab, :
名前の数が誤りです
-- &new{2006-01-03 (火) 15:41:29};
-鍋さんの悪い癖は，コメントに反応しないことだね。 -- &new{2006-01-05 (木) 01:09:07};
**R. Tierney 氏の R バイトコンパイラーについて
>[[ma]] (2005-12-29 (木) 18:12:39)~
~
何でも掲示板に出ている「R. Tierney 氏の R バイトコンパイラー（http://www.okada.jp.org/RWiki/index.php#content_1_50）」について質問させて下さい。~
~
”bad function”というエラーが出た場合には具体的にどのような問題が考えられるのでしょうか？ 附属のmanには特に関連する記述を見つけられませんでした。~
漠然とした質問で申し訳ないのですが、どなたかお教え頂けないでしょうか？~
//
-少々説明を追加させて下さい。cmpfun()でバイトコンパイルはできているようなのですが、その関数を実行すると"bad function"というエラーがでます。 -- [[ma]] &new{2005-12-30 (金) 00:57:21};
-もしあまり長くなければ、エラーが起きる実例を紹介するとコメントのしようがあるかもしれません。いずれにしても、制限はかなりあると思います。 -- &new{2006-01-01 (日) 01:31:21};
**軸の省略
>[[輝]] (2005-12-23 (金) 16:33:08)~
~
barplot()で，図を描いており，一つだけ大きい項目があり，~
途中のY軸(目盛り）を省きたいのですが，~
どのようにすればよいのですか？~
//
-[[Ｑ＆Ａ (初級者コース)#content_1_29]]にある，「ヒストグラムのy軸に中断線を入れたいです」というのと似たようなことなんでしょうか。。。。否定的。 -- &new{2005-12-23 (金) 16:58:32};
-先の質問をしたAkiraです。外れ値がある場合は対数変換という手もあると思います。それと外れ値の除外は少し難しい問題ですね（技術的という意味ではなく）。 -- [[Akira]] &new{2005-12-27 (火) 08:26:19};
**RExcelの起動時エラー
>[[YK]] (2005-12-23 (金) 12:52:29)~
~
RExcelを利用しようとしたのですが，Excel起動時にエラーが発生しました．~
~
エラーメッセージは，以下のとおりです．
「実行時エラー 5 プロシージャの呼び出し、または引数が不正です。」
そのエラーの発生箇所は，☆の付いている箇所です．
Private Sub RemoveRMenu()
On Error Resume Next
☆CommandBars(1).Controls("RExcel").Delete
End Sub
RExcelの実行した環境は以下のとおりです．
Windows XP(SP2)
Office Excel 2003
R 2.2.1
RServ 2.0.0（Rserv1.3.5でも同様にエラーが発生しました）
RExcelに添付されていたExcel 02 - Macro Demoを実行した際にも同様のエラーが発生しました．~
~
エラー修正，回避などのヒントをご教授願えないでしょうか。~
よろしくお願いいたします。~
//
-思いつきです。3年ほど前の話ですが、RExcelはエラーが多発しました。そのときはバージョンとかで愛称があるのかなと思い、それ以来使っていません。WinXP(SP1)&Office2000でした。 -- [[Akira]] &new{2005-12-23 (金) 22:21:44};
-ありがとうございます．バージョンによるエラーを視野に入れて，試していきたいと思います． -- [[YK]] &new{2005-12-26 (月) 16:23:05};
**MASS lda 判別関数と偏F値について
>[[mt]] (2005-12-21 (水) 11:44:53)~
~
MASS ライブラリ中の lda() 関数を実行した結果から、判別関数を求める方法とそれぞれの変数についての偏F値の求め方がわかりません。~
R中ではpredictで判別してくれるのでいいのですが、実際どの変数が分類に寄与していて、どういったルールで判別しているのかがわかりません。~
どなたかヒントなどご教授願えないでしょうか。~
よろしくお願いいたします。~
//
-lda は，正準判別分析を行っていますね。ということで，偏F値と言う概念はないのではないかな。あれは，フィッシャーの線形判別分析でのことですから。結果のわかっているデータを分析してみて，比べてみると良いでしょう。また，ソースを見ることができます。MASS:::lda.default などや，MASS:::predict.lda など。プログラムのソースだけを追跡するのは難しい作業なので，プログラムを自前の関数にコピーして，print関数を挿入してプログラム中の変数名（オブジェクト）が何を表しているのかを確認していくとよいでしょうね。どういう風に計算しているのかわかったら教えてください。 -- &new{2005-12-21 (水) 12:48:44};
**R2.2.0でのeps出力について。
>[[KJ]] (2005-12-20 (火) 21:37:18)~
~
Windows 版の R-2.2.0 で eps 形式での出力を行ったファイルを GSview で表示しようとすると、GSview から以下のようなエラーメッセージが吐き出されます。~
それ以前のバージョンの R であれば問題なく表示できるのですが…。~
原因は何かわかりますでしょうか？
GSview 4.7 2005-03-26
Unknown in Comments section at line 2:
%%DocumentNeededResources: font Times-Roman
Unknown in Comments section at line 3:
%%+ font Times-Bold
Unknown in Comments section at line 4:
%%+ font Times-Italic
中略
End offending input ---
file offset = 5390
gsapi_run_string_continue returns -101
~
//
-うーん, gsの環境のような気がするんですが...もし差し支えなければ問題のファイルを添付してみては? -- [[なかま]] &new{2005-12-21 (水) 13:56:50};
-ありがとうございます。
x <- seq(-3,3,by=0.01)
plot(x,dnorm(x),type="l")
として、右クリックでeps保存したファイルです。~ -- [[KJ]]
#ref(http://www.sfc.keio.ac.jp/~kj198276/test.eps)
-中は闇Ｒ-2.2.0のようですが,多分GSがCID未対応の設定もしくはバージョンなんだと思います. Rprofileに
setHook(packageEvent("grDevices", "onLoad"),
function(...) grDevices::ps.options(cidfamily=""))
等として日本語環境を無効にするか, gsが日本語を扱えるようにきちんと設定して
あげてみてください。 -- [[なかま]] &new{2005-12-21 (水) 22:51:02};
-あ、闇版を使わないと言う手もあります。(でも、2.3.0になるとやっぱりきちんと設定する必要が出てきますね) かくゆう私も、GSの設定は面倒なので...(^-^; -- [[なかま]] &new{2005-12-21 (水) 22:59:22};
-知識不足で、仰っている意味が捉えきれません…ちょっと調べてみます。 -- [[KJ]] &new{2005-12-22 (木) 01:56:34};
-やはりGSの設定の問題でした。情報ありがとうございました！ -- [[KJ]] &new{2005-12-23 (金) 12:48:30};
**barplot( )の凡例位置の変更
>[[兄歯]] (2005-12-17 (土) 14:17:53)~
~
barplot( )の凡例位置がグラフと重なってしまうので，凡例の位置を変更したのですが。~
//
-legend()はいかがでしょうか？~
x <- matrix(rnorm(100), ncol=10)
dimnames(x) <- list(1:10, 1:10)
barplot(x, beside=TRUE, legend.text=rownames(x), col=1:10)#こういうことでしょうか？
barplot(x, beside=TRUE, col=1:10)
legend("topright", legend=rownames(x), pch=15, col=1:10)#これは上と同じ
legend("topleft", legend=rownames(x), pch=15, col=1:10)#これはどうでしょうか？
-- [[Akira]] &new{2005-12-17 (土) 15:33:28};
-解決しました。どうもありがとうございます。 -- [[兄歯]] &new{2005-12-17 (土) 16:35:44};
**軸の分割数について
>[[内藤]] (2005-12-16 (金) 16:57:40)~
~
plot( )で，軸（目盛り)の分割を任意に与えたいのですが，~
どのようにすればよいのですか？~
//
-? axis をやってみればよろしいかと。
hist(rnorm(1000),xaxt="n")
axis(1,pos=0,at=c(-1,0.5,2.8),labels=TRUE)
などで，遊んでみる。 -- &new{2005-12-16 (金) 18:00:36};
-どうもありがとうございます。助かりました。 -- [[内藤]] &new{2005-12-16 (金) 20:05:16};
**maptoolsの座標の参照方法
>[[traveller]] (2005-12-06 (火) 12:42:34)~
~
maptoolsの各図形を参照する方法、どなたか、教えていただけませんか？~
//
-各図形を参照するとは？もう少し具体的に。何をしたいのかも合わせて説明してくれると答えやすい。 -- [[谷村]] &new{2005-12-06 (火) 14:12:53};
-ポイント、ラインおよびポリゴンの座標、BBOXを参照する方法をご教示くださればありがたいのですが。これらについては、maptoolsのヘルプには見当たらないようですが． -- [[traveller]] &new{2005-12-06 (火) 14:14:10};
-maptoolsで扱うObjectはMap Obejctと呼ばれ、"Shapes"と"att.data"からなります。Shapeの方にお望みの情報が格納されています。リストを理解できていれば簡単に取り出すことができます。 -- [[谷村]] &new{2005-12-06 (火) 14:20:51};
> x <- read.shape(system.file("shapes/sids.shp", package = "maptools")[1])
> x$Shapes
[[1]]
[[1]]$Pstart
[1] 0
[[1]]$verts
[,1] [,2]
[1,] -81.47276 36.23436
[2,] -81.54084 36.27251
[3,] -81.56198 36.27359
[略]
[26,] -81.45289 36.23959
[27,] -81.47276 36.23436
[[1]]$shp.type
[1] 5
[[1]]$nVerts
[1] 27
[[1]]$nParts
[1] 1
[[1]]$bbox
[1] -81.74107 36.23436 -81.23989 36.58965
[[1]]$shpID
[1] 0
attr(,"nVerts")
[1] 27
attr(,"nParts")
[1] 1
attr(,"shp.type")
[1] 5
attr(,"bbox")
[1] -81.74107 36.23436 -81.23989 36.58965
[以下略]
従って、例えば、1つ目のポリゴンの座標を得る場合には
> x$Shapes[[1]]$verts
でOKです。BBOX座標が必要な場合は、
> x$Shapes[[1]]$bbox
ですね。
-谷村さん、ご親切なご回答ありがとうございます。 -- [[traveller]] &new{2005-12-06 (火) 14:40:52};
**ベイズクラスタ
>[[ド級初心者]] (2005-11-07 (月) 12:58:10)~
~
Rでベイズクラスタ分析が可能と風の噂で聞きましたが、どのプログラム(zipファイ名)が該当するのでしょうか？また、その使用方が記述されたマニュアルに相当するものは存在するのでしょうか？あれば、サイトアドレスなどご教示下さい。~
//
**コーフェン相関係数の求め方
>[[mori]] (2005-11-02 (水) 14:08:03)~
~
階層的クラスタリングでできた樹形図（デンドログラム）の結合距離と元データの距離の相関を求めるコーフェン相関係数を求める方法はないでしょうか？または階層的クラスタリングの結果からコーフェン行列を求める方法でもいいのですが。アドバイスお願いします。~
//
-知らないものだったのでググッてみた。~
Sneath and Sokal『数理分類学』（1994年，内田老鶴圃，ISBN: 4-7536-0117-X）はお読みか？~
ひどい訳本みたいで，三中先生に「だいたいこの訳本，トンデモな新造訳語のオンパレードで 中略 共表形相関係数だって「コーフェン相関係数」（p.321）などと訳されていて，もう話にならん．〈Cophenetic〉がどーして「コーフェン」になるの？？」とか書かれている。~
[[http://cse.niaes.affrc.go.jp/minaka/diary2003-09.html:http://cse.niaes.affrc.go.jp/minaka/diary2003-09.html]]~
[[http://www.i-juse.co.jp/statistics/support/sympo/m10/m10-s01.html:http://www.i-juse.co.jp/statistics/support/sympo/m10/m10-s01.html]]に解説があったりするようだし，~
[[http://www.sci.kagoshima-u.ac.jp/~dllsa/jp/dlls/cluscrt/cluscrt.html:http://www.sci.kagoshima-u.ac.jp/~dllsa/jp/dlls/cluscrt/cluscrt.html]]には dll なんぞがあるようだ。 -- &new{2005-11-02 (水) 17:56:39};
-ちゃんとした名前がわかれば，さらにググッてみることができる。cophenetic という関数が stats にあるではないか。目的のことをやってくれるのかどうか。。。さらに調べよう。 -- &new{2005-11-02 (水) 18:19:26};
-ありがとうございます。その関数のようです。綴りはcoffenだと思ってググッてました・・・ -- [[mori]] &new{2005-11-04 (金) 08:09:43};
-"どーして「コーフェン」になるの" -> これは multicolinearity を「マルチコ」と略す類でしょうか。まさか「コーフェンさんの」と勘違いしたというわけでは。 -- &new{2005-11-04 (金) 12:32:55};
**ロジスティック回帰での変数指定方法とresidualsなどの値
>[[たかくら]] (2005-11-01 (火) 19:00:18)~
~
glmでロジスティック回帰をしていてほとんど同じ（？）計算をしているつもりなのに、residualsの値が大きく変わることに気が付きました。どういう理屈なのでしょうか？~
例えば、独立変数xが１，２，３のそれぞれの場合に、50試行のうち成功が5, 20, 45例だったとします。~
このとき、試行数や成功例数の与え方としてweightsパラメータで渡すのも、データフレームを直接放り込むのも同じだと思っていました。~
~
resp<-rep(c(1,0),3)~
x<-(1,1,2,2,3,3)~
weight<-c(5,45,20,30,45,5)~
df1<-data.frame(resp,x,weight)~
summary(glm(resp~x,weights=weight,family=binomial,data=df1))~
~
x<-(1:3)~
y<-c(5,20,45)~
df2<-data.frame(x,y)~
df2$Ymat<-cbind(df2$y,50-df2$y)~
summary(glm(Ymat~x,family=binomial,data=df2))~
~
このどちらの場合もデータとしては同じだと思いますし、推定されるcoefficientも同じようです。~
でも、residualsの値やAICなど、モデル全体の適合性に関わるパラメータは全然違います。どうしてなんでしょうか？~
//
-なぜ同じと思ったのか不思議ですね。後者は，データは3対ですよ。degrees of freedom の違いに注意。 -- &new{2005-11-01 (火) 20:10:47};
-それと，weights は本来の意味とも違うのでは？(resp=1,x=1) というデータが 5 対などなどというのとは違うでしょう。これこそ 自由度により一層大きな違いが出る。。 -- &new{2005-11-01 (火) 20:25:25};
-あなたのやりたかったことは，むしろこちらのほうなんではないか？~
> x <- rep(1:3, each=50)
> y <- c(rep(1,5), rep(0,45), rep(1,20), rep(0,30), rep(1,45), rep(0,5))
> summary(glm(y ~ x, family=binomial))
Call:
glm(formula = y ~ x, family = binomial)
Deviance Residuals:
Min 1Q Median 3Q Max
-2.0559 -0.4062 -0.4062 0.5075 2.2521
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -4.6727 0.7323 -6.381 1.76e-10 ***
x 2.2191 0.3366 6.594 4.29e-11 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 207.28 on 149 degrees of freedom
Residual deviance: 133.16 on 148 degrees of freedom
AIC: 137.16
Number of Fisher Scoring iterations: 4
違うか？ -- &new{2005-11-01 (火) 21:12:22};
-どうもありがとうございます。ご指摘の通りです。それをやりたかったんだけど、データファイルがでかくなるのが嫌で…weightsパラメータで何とかできないものかと思ったんです。
weightsは重み付けだけで、自由度が変化してないのですね。たしかに。
重み付けだけというのが具体例としてぴんと来ません。わたしが書いた２種の解析方法はどういう実験デザインに対応していると理解したらいいのでしょうか？ -- [[たかくら]] &new{2005-11-02 (水) 11:44:56};
-いくつか例を当たってみましたが、"Introductory Statistics With R"を見る限り、responceとしてベクトルを与える最初の方法も、生データを与える最後の方法もどちらもありみたいですね。weightsパラメータで重み付けを行う方法もありですが、
rate <- c(5/50, 20/50, 45/50)
total <- c(50,50,50)
x <- 1:3
summary(glm(rate~x,weights=total,binomial))
とする方が、よりスマートですね。１番目、２番目の方法では、生データを直接使う方法に比べてdfもdevianceも変わりますが、Dalgaardは気にすんなと書いてました。 -- [[たかくら]] &new{2005-11-04 (金) 15:30:04};
**有効桁数
>[[初心者K]] (2005-10-28 (金) 12:11:29)~
~
コンピュータでは2進数で数字が表現されているため、~
近似的な計算によって、桁落ちや丸め誤差が生じてしまうと聞きましたが、~
それはコンピュータでの計算に一般的に言えることなのでしょうか？~
（やはりRでもアルゴリズムによって、同様の事が生じてしまうのでしょうか。）~
Rではどのような扱いになっているのか教えてください。~
//
-整数は正確に保持されています（どんな大きな整数でもというわけではないです）。~
小数部分を持つ数も，小数部が2進数で表されるもの（0.5, 0.25, 0.125 などなどや，それらの和で表現されるもの）は，正確に保持されます。~
それ以外の数は，"近似値"が保持されますが，実際上は有効桁が22桁ほどありますので，実用上は問題ないです。~
ただ，アルゴリズムに気をつけないととんでもない結果が得られることもあります。~
よく知られていることではつい最近までの Excel の分散などを求める関数はまずいアルゴリズムに従っていたので特殊なデータの場合には不都合な結果が得られていました。~
また別の例では，二次方程式の解を求める際に，単に解の公式から二つ求めるのではなく，適切な方の一つの解は解の公式から求め，もう一つは求めた解から解と係数の関係を用いて導き出すべきだとかいうことがあります。~
よくうっかり落ち込むのは，演算結果が特定の定数値に等しいかどうかの判定を行う場合に，決して等しいかどうかを判定してはならないというようなことですね。~
このようなことは R でもその他のソフトでも，あるいは C 言語などでプログラムを書くときにも同じような問題が存在します。だって，R も何らかの言語で書かれたソフトですからね。 -- &new{2005-10-28 (金) 13:25:34};
-桁落ちは，2 進数であろうと10進数であろうと生じます。10進数で 1234.56789 - 1234.56788 = 0.00001 と9桁の有効桁を持っていたのがあっという間に1桁しかなくなってしまうのが桁落ちです。先の，二次方程式の解を求めるときにこれが生じる場合があるということ。他に，積み残しというのがある。これも5桁しか扱えない10進数コンピュータで，12345 + 0.0001 = 12345 としかなりようがない。これが，Excel の分散の計算における"ちょんぼ"の原因。 -- &new{2005-10-28 (金) 13:38:24};
-ありがとうございました。それと関連して追加の質問なのですが、Rは何バイトで計算がされているのでしょうか？（この表現方法は適当でないかもしれませんが。）教えて頂けたら幸いです。 -- [[初心者K]] &new{2005-10-28 (金) 16:45:46};
-64bit 版もあるようですし。数値演算プロセッサを使って計算したりで，このあたりの詳しいことは，なかまさんが詳しいのでしょう。 -- &new{2005-10-28 (金) 17:49:07};
-数値演算はど素人(実際そうですのでご心配なく)ですが,だいたい信用してもよさげな値は15,6桁(と言われていますが,実際はデータ次第です)です.演算はたいていdouble型(64bit)で行われています.たいていと言うのはコンパイルオプションやCPUによって一部80bit等で行われる場合もあるからです.64bit版と言うのは整数演算に関するメリットは全く無くってやや大きな配列が扱える(それが嬉しいのもある)が現在のRでは添字がint(31bit)に制限されています.この他,コンパイル時の最適化によって(演算順序やら置き換えやら沢山...)も変わります.基本的に精度を大事(ただ,精度が上がれば良いと言うわけでも無いのが痛い)にしたオプションを選ばないとRはおかしくなります.help(all.equal)等としてみると,tolerance = .Machine$double.eps ^ 0.5 等となっています. 要素が増えれば積み残しや丸め誤差が増えるでしょうから適当に判断してやる必要がある場合もあります. とここまでは,CPU,コンパイラに関する部分です. 次はアルゴリズムに関する部分です.ある問題を解くのにA,Bと二つの...と言う問題が沢山...等といった問題があるの*だけ*はなんとなく素人ながらわかりました. と言うような話をつらつら書いていけば随分な量になるんではないかと,残念ながら私の技量で書くのは不可能なのではないかと...(^_^; -- [[なかま]] &new{2005-10-28 (金) 20:48:29};
-long doubleは処理系依存です. -- [[なかま]] &new{2005-10-28 (金) 21:49:28};
-long double というのは定義上は4倍精度で34，5桁でしたっけね。C でプログラムを書くことがなくなって，とんと疎くなってしまいました。しかし long double は処理系依存で，私が使っていた BCC とか GCC では，19桁程度を使えていたようです。- &new{2005-10-28 (金) 21:35:12};
-そうはいってもアルゴリズムを気をつけなければならないということでしょうかね。ちなみに64bitというのは何か他のものと比較したときに、どの程度のレベルなのでしょうか？ -- [[初心者K]] &new{2005-10-29 (土) 11:35:28};
-「何か他のもの」<- 意味不明で良い質問ではないですね。倍精度計算の次は Mathematica 等の任意精度計算になってしまうのでは。中間さんが書かれているように、この問題はあまりに重大過ぎて、従って「普通見てみぬふりをする」類の問題です。まず数値解析の入門書を一冊読み、そして自分で書いたプログラムのおかしな結果に悩みまくって、勘どころ、妥協点を体得する他無いでしょう（と経験から思います）。64 ビットで十分とは今や誰も思っていないでしょうが、これまでの（64 bit用）ソフトの膨大な蓄積を考えると誰も手が付けられない、というのが真相だと思います。 -- [[QDU]] &new{2005-10-29 (土) 13:55:47};
-ありがとうございました。頂いたコメントを参考に、またさらに悩んで、体得したいと思います。 -- [[初心者K]] &new{2005-11-07 (月) 00:05:45};
**ARMAモデルの最尤推定
>[[tanabe]] (2005-10-27 (木) 23:30:27)~
~
いつもお世話になっております。~
~
Rでは時系列のARおよびARIMAモデルに関してパラメータの最尤推定を行うことができますが、このアルゴリズムはどのようなものなのでしょうか？正確なのか擬似なのか等知っておられる方がいらっっしゃいましたら教えてください。~
また、多変量ARMAモデルのパラメータ最尤推定に対応した機能はありますか？~
さらに、欠損値を含むデータに対しても対応する機能があるのか知っておられるようでしたら教えてください。~
~
以上、よろしく御願いします。~
//
**時系列の成分分解
>[[kkondo]] (2005-09-29 (木) 13:33:10)~
~
１変量時系列データを、2階型トレンド成分、ARMA(p,q)成分、白色雑音成分の和へ分解する計算を行ないたいと思っておりますが、R の decompose, stl 関数はこのような分解能力はなさそう(?)です。~
最終的には、AIC (もしくはBIC) が最適になるように、成分分解、すなわち、各成分の振幅を与える分散の大きさと、ARMA(p,q) の次数と係数の値を計算したいと思っておりますが、どこかに参考情報はあるでしょうか？~
//
-自己補足説明なのですが、標準的(?)な方法論を用いた方法を念頭においております。たとえば、G. Kitagawa, T. Higuchi, and F. N. Kondo, Smoothness Prior Approach to Explore Mean Structure in Large-scale Time Series, Theoretical Computer Science 292, No.2, 431-446, 2003. http://tswww.ism.ac.jp/higuchi/index_e/papers/KHKforTCSfin.pdf -- [[kkondo]] &new{2005-09-29 (木) 13:59:56};
**プロットで使える文字
>近藤 (2005-09-28 (水) 23:47:54)~
~
expression() で使える文字の一覧というのは存在するのでしょうか？~
expression(cdots) は存在するのですが expression(cdot) はないんですね。~
//
-help("plotmath")に一覧があります.demo(plotmath)で一覧が描画されます. -- [[なかま]] &new{2005-09-29 (木) 00:32:02};
-印刷される場合は[[R-Tips:http://cse.naro.affrc.go.jp/takezawa/r-tips.pdf]]の133ページ〜138ページをごらんください． -- &new{2005-09-29 (木) 00:44:39};
-ありがとうございました。 中点というのは存在しないみたいですね。 tex の $x ?cdot y$ です。 -- 近藤 &new{2005-09-29 (木) 01:36:46};
**R-2.1.1 (Windows版)のコマンドラインからインストール
>[[BR103]] (2005-09-22 (木) 16:28:39)~
~
題名のように、GUIを使わずにインストールはできないのでしょうか？~
//
-rw2011.exe /VERYSILENT しかし、これはデフォルトでインストールされる. src/gnuwin32/installer/ のファイルをInnoのhelpと見比べてデフォルトで目的の環境を構築出来るようにするのが良いと思う.インストーラだけの再構築なら,InnoとRのソースとActivePerlだけで可能なはず.(手間を減らすには手間をかける) -- [[なかま]] &new{2005-09-23 (金) 01:12:54};
**read.dbfの不具合？
>[[谷村]] (2005-08-19 (金) 16:11:14)~
~
r-helpでスルーされたので、こちらに相談します。~
KEYCODE (数値型19桁小数0桁)
422010010
42201002101
42201002102
42201002103
42201002104
422010060
422010071
422010072
42201008001
42201008002
という内容のtest.dbfとtest.csvを用意します。&size(8){実物を添付しようとしたのですが、画像ファイルでないということで添付できませんでした。};~
> library(foreign)
> cbind(read.csv("test.csv"),read.dbf("test.dbf"))
KEYCODE KEYCODE
1 422010010 422010010
2 42201002101 NA
3 42201002102 NA
4 42201002103 NA
5 42201002104 NA
6 422010060 422010060
7 422010071 422010071
8 422010072 422010072
9 42201008001 NA
10 42201008002 NA
となります。read.shapeで読み込んだ属性値がおかしいので、調べたら read.dbfが期待通りの結果になっていないことに気がつきました。~
OpenOffice.org Calcまたはtxt2dbfでtest.dbfを作成して試しました。どちらも同じ結果です。こちらの環境は、R-2.1.1+foreign 0.8-9です。どなたか追試験だけでも、お願いします。~
//
-dbfファイルとは、dBASEの形式でしょうか？ -- &new{2005-08-19 (金) 16:16:52};
-dbfはdBASE形式ファイルの拡張子です。 -- [[谷村]] &new{2005-08-19 (金) 16:20:00};
-write.dbf で書き出したものとは違うのですね。そちらが作った test.dbf なりをアップロードしておいてくれれば，追試可能だとは思いますが。 -- &new{2005-08-19 (金) 16:40:18};
-添付できないので、私が管理するサーバにアップロードしました。&ref(http://shakan2.tm.nagasaki-u.ac.jp/~umusus/test.csv); &ref(http://shakan2.tm.nagasaki-u.ac.jp/~umusus/test.dbf); 明日(8/20)は停電のためサーバがダウンします。 -- [[谷村]] &new{2005-08-19 (金) 16:46:19};
-セキュリティホール対策のため、画像以外の添付機能は無効にしてあるのでご不便をおかけして申し訳ありません。R(MacOS X) 2.1.1 + foreign 0.8-8 ですが、
> cbind(read.csv("t1.csv",header=F),read.dbf("t2.dbf"))
V1 N1
1 422010010 422010010
2 42201002101 42201002101
3 42201002102 42201002102
4 42201002103 42201002103
5 42201002104 42201002104
6 422010060 422010060
7 422010071 422010071
8 422010072 422010072
9 42201008001 42201008001
10 42201008002 42201008002
と問題ないようです。dbfはNeoOffice/J 1.1 RC Patch8で作成しましたが、NeoOfficeの動作はほぼOpenOffice.orgと同様と思います。 -- [[岡田]] &new{2005-08-19 (金) 16:45:48};
-上記からダウンロードした test.dbf を読み込んだところ、最初の出力結果でNAのところが手元の環境でも 2147483647 となって異常です。dbfファイルがおかしいのかもしれません。 -- [[岡田]] &new{2005-08-19 (金) 16:52:02};
-予想通り。4バイト整数として読んでいるようなので，Mac では以下のようになります。~
> cbind(read.csv("test.csv"),read.dbf("test.dbf"))
KEYCODE KEYCODE
1 422010010.00000000 422010010
2 42201002101.00000000 2147483647
3 42201002102.00000000 2147483647
4 42201002103.00000000 2147483647
5 42201002104.00000000 2147483647
6 422010060.00000000 422010060
7 422010071.00000000 422010071
8 422010072.00000000 422010072
9 42201008001.00000000 2147483647
10 42201008002.00000000 2147483647
> x <- cbind(read.csv("test.csv"),read.dbf("test.dbf"))
> class(x[2,2])
[1] "integer"
岡田さんと同じ環境なんだけど，ファイルは，谷村さんのをダウンロードして使いました。 -- [[青木繁伸]] &new{2005-08-19 (金) 16:50:04};
-ご無沙汰です。> 岡田さん 追試験をありがとうございました。こちらでは、txt2dbfでもoocalcでもArcGISでも再現します。dBASEのバージョンの問題でしょうか。
$ dbf2txt test.dbf
#KEYCODE
422010010
42201002101
42201002102
42201002103
42201002104
422010060
422010071
422010072
42201008001
42201008002
dbfファイルは正常に見えます。oocalcで開いても正常に開きます。 -- [[谷村]] &new{2005-08-19 (金) 16:52:55};
-signed integer の桁あふれみたいですね。2147483647 = 0x7FFFFFFF -- [[岡田]] &new{2005-08-19 (金) 17:01:17};
-ちなみにppc64 Rでも 2147483647 でした... foreignのソースコードをざっとみると32ビット決めうちのような感じです。 -- [[岡田]] &new{2005-08-19 (金) 17:07:14};
-正常に読めたdbfは [[ここ:http://epidemiology.md.tsukuba.ac.jp/~mokada/R/t2.dbf]]です。 -- [[岡田]] &new{2005-08-19 (金) 17:09:36};
-岡田さんも青木先生もありがとうございます。原因が絞り込まれて大変助かります。実際に扱っているdbfはshpファイルの一部なので、dbfに手を加えるのはちょっと厳しいです。何かR側で打つ手はありますでしょうか。 -- [[谷村]] &new{2005-08-19 (金) 17:15:16};
-下記のような試行をしてみました。-- [[谷村]] &new{2005-08-19 (金) 17:16:06};
> write.dbf(read.csv("test.csv"),file="test2.dbf")
> cbind(read.csv("test.csv"),read.dbf("test2.dbf"))
KEYCODE KEYCODE1
422010010 4220100102
42201002101 422010021013
42201002102 422010021024
42201002103 422010021035
42201002104 422010021046
422010060 4220100607
422010071 4220100718
422010072 4220100729
42201008001 42201008001
42201008002 42201008002
-岡田さんのと谷村さんのでは，同じ dBase 形式でも中身が違いますね。岡田さんの方のファイルは小数点が入っているし。。~
$ > hexdump t2.dbf
0000000 0305 0813 0a00 0000 4100 0f00 0000 0000
0000010 0000 0000 0000 0000 0000 0000 0000 0000
0000020 4e31 0000 0000 0000 0000 004e 0000 0000
0000030 0e02 0000 0000 0000 0000 0000 0000 0000
0000040 0d20 3432 3230 3130 3031 302e 3030 0000
0000050 2034 3232 3031 3030 3231 3031 2e30 3020
0000060 3432 3230 3130 3032 3130 322e 3030 2034
0000070 3232 3031 3030 3231 3033 2e30 3020 3432
0000080 3230 3130 3032 3130 342e 3030 2034 3232
0000090 3031 3030 3630 2e30 3000 0020 3432 3230
00000a0 3130 3037 312e 3030 0000 2034 3232 3031
00000b0 3030 3732 2e30 3000 0020 3432 3230 3130
00000c0 3038 3030 312e 3030 2034 3232 3031 3030
00000d0 3830 3032 2e30 3000
00000d7
$ > hexdump test.dbf
0000000 0369 0813 0a00 0000 4100 1400 0000 0000
0000010 0000 0000 0000 0000 0000 0000 0000 0000
0000020 4b45 5943 4f44 4500 0000 004e 0000 0000
0000030 1300 0000 0000 0000 0000 0000 0000 0000
0000040 0d20 2020 2020 2020 2020 2020 3432 3230
0000050 3130 3031 3020 2020 2020 2020 2020 3432
0000060 3230 3130 3032 3130 3120 2020 2020 2020
0000070 2020 3432 3230 3130 3032 3130 3220 2020
0000080 2020 2020 2020 3432 3230 3130 3032 3130
0000090 3320 2020 2020 2020 2020 3432 3230 3130
00000a0 3032 3130 3420 2020 2020 2020 2020 2020
00000b0 3432 3230 3130 3036 3020 2020 2020 2020
00000c0 2020 2020 3432 3230 3130 3037 3120 2020
00000d0 2020 2020 2020 2020 3432 3230 3130 3037
00000e0 3220 2020 2020 2020 2020 3432 3230 3130
00000f0 3038 3030 3120 2020 2020 2020 2020 3432
0000100 3230 3130 3038 3030 3200
0000109
dBase ファイルは比較的おとなしい構造をしている（上記）ので，むかし，C で読み取るプログラムを書いたことがある。 -- [[青木繁伸]] &new{2005-08-19 (金) 17:21:48};
-単純に scan で読んで，その後選択変換する。。。~
> x <- scan("test.dbf", what="")
> x
[1] "?003i?b?023" "" "A" "422010010"
[5] "42201002101" "42201002102" "42201002103" "42201002104"
[9] "422010060" "422010071" "422010072" "42201008001"
[13] "42201008002"
> as.numeric(x[-(1:3)])
[1] 422010010.00000000 42201002101.00000000
[3] 42201002102.00000000 42201002103.00000000
[5] 42201002104.00000000 422010060.00000000
[7] 422010071.00000000 422010072.00000000
[9] 42201008001.00000000 42201008002.00000000
などとねぇ。 -- [[青木繁伸]] &new{2005-08-19 (金) 17:26:58};
-数値型19桁小数0桁ではなく数値型19桁小数5桁にしてみると、桁あふれが無くなりました。でも、ArcGISが吐くシェープファイルに含まれるdbfを読み込みたいので、何とかdbfをさわらずにできないものか。。。。 -- [[谷村]] &new{2005-08-19 (金) 17:34:10};
-なるほど、scanを使う方法があるのですね。read.shapeの改造にチャレンジしてみます。 -- [[谷村]] &new{2005-08-19 (金) 17:37:34};
-無理矢理 foreign のソース自体を書き換えるという方法もあります。ちょっと強引ですが(^^; -- [[岡田]] &new{2005-08-19 (金) 17:41:05};
-うんと,こんなんでいいんじゃないかなぁ -- [[なかま]] &new{2005-08-19 (金) 17:45:12};
--- foreign.orig/src/dbfopen.c 2004-10-27 20:26:43.000000000 +0900
+++ foreign/src/dbfopen.c 2005-08-19 17:42:47.549620408 +0900
@@ -970,7 +970,8 @@
|| psDBF->pachFieldType[iField] == 'F' )
/* || psDBF->pachFieldType[iField] == 'D' ) D is Date */
{
- if( psDBF->panFieldDecimals[iField] > 0 )
+ if( psDBF->panFieldDecimals[iField] > 0 ||
+ psDBF->panFieldSize[iField] > 9 )
return( FTDouble );
else
return( FTInteger );
- かぶった(^^;
*** foreign/src/Rdbfread.c 2005-04-06 04:35:39.000000000 +0900
--- foreign.new/src/Rdbfread.c 2005-08-19 17:39:25.000000000 +0900
***************
*** 89,95 ****
SET_VECTOR_ELT(df, nRvar, allocVector(STRSXP,nrecs));
break;
case 2:
! SET_VECTOR_ELT(df, nRvar, allocVector(INTSXP,nrecs));
break;
case 3:
SET_VECTOR_ELT(df, nRvar, allocVector(REALSXP,nrecs));
--- 89,95 ----
SET_VECTOR_ELT(df, nRvar, allocVector(STRSXP,nrecs));
break;
case 2:
! SET_VECTOR_ELT(df, nRvar, allocVector(REALSXP,nrecs));
break;
case 3:
SET_VECTOR_ELT(df, nRvar, allocVector(REALSXP,nrecs));
***************
*** 120,129 ****
case 2:
if( DBFIsAttributeNULL( hDBF, iRecord, i ))
! INTEGER(VECTOR_ELT(df, nRvar))[iRecord] = NA_INTEGER;
else
! INTEGER(VECTOR_ELT(df, nRvar))[iRecord] =
! DBFReadIntegerAttribute( hDBF, iRecord, i );
nRvar++;
break;
--- 120,129 ----
case 2:
if( DBFIsAttributeNULL( hDBF, iRecord, i ))
! REAL(VECTOR_ELT(df, nRvar))[iRecord] = NA_INTEGER;
else
! REAL(VECTOR_ELT(df, nRvar))[iRecord] =
! DBFReadDoubleAttribute( hDBF, iRecord, i );
nRvar++;
break;
-- [[岡田]] &new{2005-08-19 (金) 17:49:05};
-dBaseの仕様だと要するに10進数で,Rの方は整数型と浮動小数点型しか無いんで.更にこれだと丸められてひどい目にあうかもしれませんが... -- [[なかま]] &new{2005-08-19 (金) 17:49:27};
-う、輻輳。さー、どっちがいいでしょう。(笑) -- [[なかま]] &new{2005-08-19 (金) 17:54:35};
-とりあえず、なかまさんの方を試しました。うまくいきました。冗長なのでここにペースとしませんが、感謝感謝です。パッチがうまくあたらなかったので手で修正を加えました。authorsへのフィードバックはどうしましょう。 -- [[谷村]] &new{2005-08-19 (金) 18:18:31};
-ごった煮パッケージですから,普通にSummaryとして投稿しとけば, 中の人が拾ってくれると思います. -- [[なかま]] &new{2005-08-19 (金) 18:43:23};
-岡田さんのパッチも試しました。こっちは素直にパッチがあたってくれました。こちらもうまくいきました。感謝感謝です。 -- [[谷村]] &new{2005-08-19 (金) 18:44:54};
-なかまさんのパッチのほうがスマートですから、そちらをフィードバックされるのがよいと思います。 -- [[岡田]] &new{2005-08-19 (金) 18:45:08};
-r-helpにサマリを投稿しました。ありがとうございました。 -- [[谷村]] &new{2005-08-20 (土) 02:13:55};
**A I D S モデル
>[[sonkun]] (2005-08-09 (火) 17:29:30)~
~
Rで需要分析のAlmost Ideal Demand System（A I D S モデル）を試されている方、いませんでしょうか。SASのプログラムhttp://support.sas.com/rnd/app/examples/ets/elasticity/をRに直す形で取り組んでいるのですが、SASでの回帰分析"proc syslin itsur…"のところからが上手く置き換えられません。~
どなたか詳しい方いらっしゃいましたらご助言よろしくお願いいたします。
//
-トップページにリンクのある R site search を使いキーワード "A I D S" で検索されてみたら。 -- [[QDU]] &new{2005-08-09 (火) 19:56:21};
-パッケージmicEconに A I D S 関連の関数があります。 -- &new{2005-08-18 (木) 16:32:54};
**GARCHのパラメータ推定について
>[[松本]] (2005-07-23 (土) 15:19:13)~
~
以前も時系列分析について質問させて頂きました。~
GARCHモデルを使って、時系列分析をしようとしているのですが、エラーが出て出来ませんでした。具体的には、fSeriesをダウンロードして、garchFitを使いました。エラーの内容ですが、文字化けしていてわかりませんでした。~
Rのバージョンは、2.1.1です。~
ただ、garchFit(データ、order=c(0,1))のように、ARCH推定は可能なのです。~
私の仕方が悪いのかさえ分かりませんので、教えていただけませんか？~
また、ARMAモデルの次数設定の方法を以前教えて頂いたのですが、それを用いてGARCHの次数設定も可能なのでしょうか？~
用いた導出法は、
~
for (j in 1:5)
for (i in 1:5) {
tmp <- arima(data, order=c(j, 0, i))
cat("order1=", j, "order2=", i, "AIC=", tmp$aic, "?n")
}
//
-実際に操作した内容を，全部！！！コンソールからコピー・ペーストしてください。「GARCHモデルを使って、時系列分析をしようとしているのですが、エラーが出て出来ませんでした。具体的には、fSeriesをダウンロードして、garchFitを使いました。」では，さっぱりわかりません。 -- [[青木繁伸]] &new{2005-07-23 (土) 19:21:58};
- [[ヘルプ]]を参照して，美しく投稿しましょう。
-貼り付けたエラーメッセージは，文字化けしていませんね。Rコンソールで文字化けしていると言うことは，Windows ですか？Windows で日本語が文字化けするときは，インストールについてのページを読んでちゃんとインストールすること。それと，文字化けしていないメッセージを見て，原因は推定できませんか？AIC= の後にtmp$aic が書かれるはずと思っているのですがそれが書かれていませんね。garchFit が失敗しているからでしょう。時系列データがNAを含んでいるよと言われていますが大丈夫ですか。ウォーニングは負の平方根を取ろうとしたよということで pred$e の実際の値をチェックしてみてください。 -- &new{2005-07-24 (日) 11:08:05};
-rate というオブジェクトは，あなたが定義したもの？わたしが，書かれたプログラムを実行してみたら，「以下にエラーinherits(x, "data.frame") : オブジェクト "rate" は存在しません」と怒られちゃいました。 -- &new{2005-07-24 (日) 11:10:05};
-そうです。私が使っているのは、Windowsです。~
文字化けについては、ご指摘のように対処してみます。~
ただ、文字化けしていないメッセージを見ても分かりません。tmp$aicが書かれていないのは、
気づきました。以前、ARIMAで同様の事をした際は、きちんとできたのですが。~
どうも扱うデータがおかしいのかもしれないので、もう一度チェックしてみます。~
しかし、なぜGARCH推定はだめで、ARCH推定は可能であったのでしょうか？
-rateというのは確かに私が定義したものです。具体的には、FRBに掲載されていた日次の円ドル為替レートを差分操作で直したものです。-- [[松本]] &new{2005-07-24 (日) 12:10:40};
-何でかは，私にはわかりません。~
いずれにせよ，あなたには分かるべきですが，rate というものがどういうものか分からないと，ほかの Windows R 使いも，コメントのしようがなくて，じりじりしているのでは？~
一応，エラーの原因と思われるオブジェクトが示されたのですが，そのオブジェクトがなんでそんな，変な値を持っているのか調べるのが先決では？？？？？？？？？？（幾つクエスチョンマークを付ければいいですか） -- &new{2005-07-24 (日) 19:15:25};
-rateというのは、ｔｘｔ形式のファイルです。 -- [[松本]] &new{2005-07-24 (日) 20:23:35};
-t x t 形式のファイルったって，それがどのような情報を持つとお思いですか。~
極端なことを言えば，ゼロ行列かもしれないし，単位行列かもしれないわけですよ。無意味な乱数行列であってもいいのですか。そんなはずなでしょう。あなたが適切なオブジェクトを設定しているかどうか，誰も分からないのです。~
先にもコメントしたように，誰も，あなたと同じ条件で追試できないので，どこが間違いと指摘できないのです。答えがいらないのならそれで良いのですが。~
rate を作った手順そのものも，R コンソールからコピー・ペーストすればいいのです。~
要するに，あなた以外の人がその人の環境下で，あなたと同じデータを使って，同じ操作をして，同じエラーが出たときに，その人がどこが悪いかを調べて，あ，ここはこうしないとだめなんだよと指摘できるようにしないと，あなたの望む即効薬は得られないのです。~
もっとも，みんなが暇で，あなたの問題につきあってあげようと思えばのことですが。~
そうでなくても，ここが問題じゃないの？くらいのアドバイスは得られると思いますが，それも，あなたが与えられる情報次第なんですよね。たぶん。私はそうおもいますが。 -- &new{2005-07-24 (日) 20:34:08};
-答えが欲しいなら，回答があるかないか分からなくても，なるべく頻繁に回答の有無をチェックするのが礼儀でしょうね。~
レスポンス返しても，それに対する応答が長い間ないと，もうどうでも良いのかと思われませんか。 -- &new{2005-07-24 (日) 20:53:15};
-頻繁に回答のチェックはしています。 -- [[松本]] &new{2005-07-24 (日) 20:54:03};
-ただ、ここの編集機能の使い方が不慣れで、rateの説明を編集しているうちに、皆様からの返答があるもので･･･。すぐに、rateの説明を編集しますので、少し待って頂けませんか？ ~
rateですが、txtファイルに数値を打ち込んだだけのものです。
全ての数値を打ち込めないので、部分的に挙げます。
126.01
126.24
126.06
126.23
127.38
127.57
127.73
128.14
128.52
129.58
129.95
130.52
131.27
131.18
130.92
132.08
131.54
131.1
というように数値が800個程縦にならんでいます。それ以外には、行番号や、その他の数値は
一切ありません。このtxtファイルのファイル名は、rate6.txtです。
というようにしました。長く時間をかけてすみません。-- [[松本]] &new{2005-07-24 (日) 20:57:58};
-秘密事項でなければ，そのファイルを添付すればいかがですか。ファイルの添付も[[ヘルプ]]に書いてあります。 -- &new{2005-07-24 (日) 21:25:22};
-どうも
![[RjpWiki]](image/../Rlogo.png "[RjpWiki]")
