myscratchpad

基本の初歩の入門
L(p) = p³(1-p)^(5-3)
のL(p)を最大化するpを求める
# よくある例: 5回コインを投げて3回表、2回裏が出る確率 (L(p)) を最大化する、コインの表の出る確率 (p) は？

# Excelのソルバー
A2に初期値として0.5を入れる
B2に式を書く
=(A2^3)*((1-A2)^2)
ソルバーを起動し、
目的セル: $B$2
目標値: 最大化
変化させるセル: $A$2
制約条件: $A$2<=1, $A$2>=0
として実行をクリック

A2 (p) は0.5、B2 (L(p)) は0.03456が記入される


# 同じことをRでやるにはoptim関数を使う
fc <- function(x){
(x^3)*((1-x)^2) # 関数定義。上の=(A2^3)*((1-A2)^2)　と同じ
}
optim(par=0.5, fn=fc, control=list(fnscale=-1))
## parは初期値、fnは関数、controlは制御構文。fnscale=-1は最大化しろということ
## 返り値でconvergenceが0だと収束したということ

# 追記。Excelのソルバーは準ニュートン法か共役勾配 (傾斜) 法を用いている。オプションで指定できる。
# Rではmethodで指定する。大概制約条件があるので、"L-BFGS-B" であろう (デフォルトは"Nelder-Mead") 。Rjpwikiの記事を参考。感謝
optim(par=0.5, fn=fc, method="L-BFGS-B", upper=1, lower=0, control=list(fnscale=-1))


# optimize関数もある
optimize(f=fc, interval=c(0,1), maximum=T)
## fは関数, intervalは範囲、maximumは最大化

psychパッケージの出力は小数点2位までしか出力されない
psychパッケージのオブジェクトを出力するprint.psychの指定でdigits=2が強制指定されているからのようだ
head(print.psych, 5) # で見るとわかる
しかもこれはoptionsの指定を受け取らない

もっと出力したいときは別の関数でclassを変える
dat <- matrix(rnorm(100), 10, 10)
library(psych)
des.dat <- describe(dat)
class(des.dat) # "psych"クラスがつくので小数点2位までしかでない
des.dat

des.dat2 <- data.frame(describe(dat))  # data.frameでclassを変える
class(des.dat2) # "psych"クラスが外れて"data.frame"になるので、optionsに従った出力になる
des.dat2

# sample関数の覚書
sample(x=ベクトル, size=総数, replace=F:同じ要素が選ばれても良いか, prob=NULL:確率ベクトル
## 5件法のデータを100個生成する
sample(x=5, size=100, replace=T)
## prob引数。確率ベクトルとして3が半分とかのデータにする。デフォルトは一様
p <- c(0.1, 0.1, 0.5, 0.1, 0.2)
table(sample(x=5, size=100, replace=T, prob=p)) # tableで集計

# matrix関数の覚書
matrix(data=NA: データベクトル, nrow=1: 行数, ncol=1: 列数, byrow=F:
一番左の列の上から順にデータを入れていく。Tだと一番上の行の左から順にデータを入れていく, dimnames=NULL: 行列の名前をリストで与える
## 2行3列の行列をつくる。中身はNA
matrix(NA, nrow=2, ncol=3)
## 60項目の質問紙に200人が答えた、という仮想データ
## sample関数で1 から5までの数値を12000個生成し、200行60列の行列にする
mat <- matrix(data=sample(x=5, size=12000, replace=T), nrow=200, ncol=60)
summary(factor(mat)) # 集計

ずいぶん前のものみたいだけど、変なの見つけた。
何という才能の無駄使いｗ…コード欲しいな

【初音ミク】えれくとりっく・えんじぇぅをRで描いてみた



以下はリンクのみ。ニコ動のアカウント必須。
Ｒで超時空シンデレラ、ランカさんを描いてみた

全てRjpwikiの記事から。感謝。

# たくさんのグラフを一度に書くときに画面を自動的に分割する。元記事リンク
par (mfrow=n2mfrow(N)) # Nはグラフの数を入れる

# データフレームの変数名を短くする。元記事リンク
abbreviate(names(iris))
tolower(abbreviate(names(iris)))　# 小文字にするにはtolower関数を使う。大文字はtoupper

# Na, NaN, Infも集計する。元記事リンク
summary(x) で集計できる。利用例