RでHurwitz's zeta関数 フルビッツゼータ関数の近似 よいライブラリがないので、とりあえず、つくってみた。 library("Rmpfr") Hurwit_zeta_approx<-function(beta4,x){ #N<-10^7 #k<-0:N #s<-beta4 ans2<-NULL if(beta4==beta4){ for(j in 1:length(x)){ a3…

Rでパラメータ制約付非線形最適化 Rsolnpパッケージをもちいる。Rで非線形制約条件付非線形最適化を行う方法 by Rsolnp http://d.hatena.ne.jp/teramonagi/20091217/1261048574を参考に とりあえず最小二乗法でためす. library(Rsolnp)#答えと問題のデータを…

rownames(r_d)<-iconv(rownames(r_d),to="utf-8",from="cp932") colnames(r_d)<-iconv(colnames(r_d),to="utf-8",from="cp932")

Step1: ネットのつながるパソコンで wget https://cran.r-project.org/src/contrib/Archive/lavaan/lavaan_0.5-17.tar.gz 0.5-17がR 3.0で使える一番新しいものStep2: インストールしたいパソコンにコピーして以下を実行 R CMD INSTALL *0.5*17*tar.gz ERROR…

Rで因子分析と分散共分散分析をする。 因子分析は、潜在要因をみつける方法。 主成分分析と回帰分析と異なり観測変数xを説明する潜在変数をみつける方法（下のx1とx2を見つける） y1=a*x1+a*x2 y2=a*x1+a*x2 y3=a*x1+a*x2 y4=a*x1+a*x2 心理学等では、「瞬発…

dplyr ですべての列の平均を計算 dplyr::summarise_each(data,funs(mean),everything())

dplyr ですべての列の平均を計算 dplyr::summarise_each(data,funs(mean),everything())

Rでキー列ごとに代表の値をえたデータフレームをえる方法。ユニークになるように集約 つまり、グループ代表値以外データフレームから削除する方法は、distinctを用いる． http://a-habakiri.hateblo.jp/entry/2016/11/29/215013 例：入力 a b c 1 1 3 1 1 4 …

Rでカテゴリ変数の相関を調べる。質的変数とクラメール相関係数の計算。 基本は、（１）分割表を作る、（２）相関係数を調べるの流れ。 library("vcd) ##サンプルデータの作成 #x<-rbinom(1000,size=1,prob=0.5) #y<-rbinom(1000,size=1,prob=0.5+0.3*x) #x<…

ggplotでグラフの形を変える。 scale_shape_manualを用いる ggplot(smp, aes(x = year, y = value, group = country)) + geom_line(aes(color = country)) + geom_point(aes(color = country, shape = country)) + scale_shape_manual(values = 0:length(uni…

Rでエントロピーの概算： x<-rnorm(10^5) v<-hist(x,32) sum(log(v$density)*v$density*mean(diff(v$mid)),na.rm=T) 非連続なら entropy.empirical(table(x)) entropyパッケージ

"name1","name2","name3"の行だけをとりだす。 g<-c("name1","name2","name2","name3","name6","name4","name6") x<-c(3,4,5,2,1,4,5) y<-3*x+rnorm(length(x)) r1<-data.frame(x,y,g) list1<-c("name1","name2","name3") r1b<-dplyr::filter(r2,r2$g %in% l…

Rでカテゴリカルデータの相関、連関解析 https://www1.doshisha.ac.jp/~mjin/R/45/45.html http://monge.tec.fukuoka-u.ac.jp/r_analysis/data_crosstable00.html assocstats(xtabs(~y+SAIKOU,data=data,drop.unused.levels=TRUE)) コレスポンディング分析 h…

Rで２点鎖線 dashed-double-dotted line p=p+geom_vline(xintercept=3,colour=4,lwd=1,linetype="531313")linetype を文字表示でできる 531313 5の長さの線引く、３開ける、１の長さの線引く、３開ける、１の長さの線引く、３開ける 1F 1線をひく16進数のF(1…

Rで２つのカテゴリカルデータ（質的データ）の比較カイ２乗検定の例 > v1 v2 chisq.test(cbind(v1,v2">*1 Pearson's Chi-squared testdata: cbind(v1, v2) X-squared = 18.232, df = 2, p-value = 0.0001099> v2 v1 chisq.test(cbind(v1,v2">*2 Pearson's Ch…

基本は、convolve関数。線形畳み込み（普通の畳み込み）にするためには、type="o"とする。 デフォルトは循環畳み込み． q1 q2 q3 たくさん畳み込みをするときは、fft関数を利用． 0をデータの後ろにたくさん作ると循環畳み込みが線形畳み込みと同値となる。 …

Rで偏相関 library("psych") qb_learn qb_testqb_learn% dplyr::select(y,TIKUNENSUU,TINRYOU) qb_test% dplyr::select(y,TIKUNENSUU,TINRYOU)> cor(qb_learn$TINRYOU,qb_learn$y) [1] 0.1614954 > cor(qb_learn$TINRYOU,qb_learn$TIKUNENSUU) [1] -0.335243…

元データ: qba_all： y MOYORIEKI 0 大津駅 1 山科関 0 石山駅 1 大津駅 0 彦根駅 0 南草津駅 1 彦根駅 計算: 因子の順番を駅の順番を揃える stalist<-c("山科駅","大津駅","膳所駅","石山駅","瀬田駅","南草津駅","草津駅","栗東駅","守山駅","野洲駅","篠…

元データ: qba_all： y MOYORIEKI 0 大津駅 1 山科関 0 石山駅 1 大津駅 0 彦根駅 0 南草津駅 1 彦根駅 計算 qba<-qb_all #因子の順番を駅の順番を揃える #stalist<-c("山科駅","大津駅","膳所駅","石山駅","瀬田駅","南草津駅","草津駅","栗東駅","守山駅",…

Rでカテゴリカルデータでの勉強（質的データ）；数量化１類 （１）順番がつけられない質的データを考える。 それを回帰する場合適切コーディングが求められる。 ＊コーディングによって回帰係数がかわり、結果の解釈がかわる。 (A）ダメなコーディング例1: …

カテゴリカルデータで制約付の回帰。結局結論をまとめると y=a1*x1+a2*x2+a3*x3+Cで、 x=[1,2,3], x1,x2,x3ダミー変数のとき(0,1)、 a a[1] a[2] a[3] aa y x lm(y~factor(x),list(x = contr.sum)) 結果 Call: lm(formula = y ~ x, contrasts = list(x = con…

RでKLダイバージェンス（仮） 0を除かないと発散するので、除くのがポイント for(i in 1:length(list1)){ list_name<-list1[i] v1<-qb[[list_name]] a1<-v1[qb$y==1] a2<-v1[qb$y==0] aa1<-table(a1) aa2<-table(a2) ll<-unique(c(names(aa1),names(aa2))) a…

Ｒでボックスわけして最頻値 mode2<-function(x){ if(length(x)==1){ r<-x }else{ qq<-density(x,na.rm=T) qq<-approx(qq$x,qq$y,seq(quantile(qq$x,0.01),quantile(max(qq$x),1.0),length.out=1000)) r<-qq$x[max(qq$y)==qq$y] } r } #exp(log(xb,2)) xx_f<…

まずは、ポアソン過程の平均値の推定 $ cat poi.stan data{ int J; // number of data int y[J]; } parameters{ real mu; }model{ y[J] ~ poisson(mu); } 平均３０のポアソン分布 ###Poison library("rstan") y data fit iter = 1000,chains = 4) print(fit)…

x軸で等間隔． discretize関数を使う． library(infotheo) > x<-rnorm(1000) > y<-2*x+0.01*rnorm(1000) > l1<-discretize(x,disc="equalwidth",nbins=15) > tapply(x,l1[[1]],"mean") 1 2 3 4 5 6 -2.87785661 -2.50220873 -2.05626538 -1.63706018 -1.2338…

RのMCMCライブラリのRstanの導入 windows 8 64bit R x64 3.0.2 これに従う． https://github.com/stan-dev/rstan/wiki/RStan-Getting-Started （１）基本の整備 バージョンにあったRtoolsのインストール （２）gccのパスをRtoolsのものにする（３）インスト…

RのMCMCライブラリのRstanの導入 windows 8 64bit R x64 3.0.2 これに従う． https://github.com/stan-dev/rstan/wiki/RStan-Getting-Started （１）基本の整備 バージョンにあったRtoolsのインストール （２）gccのパスをRtoolsのものにする（３）インスト…

RのMCMCライブラリのRstanの導入 windows 8 64bit R x64 3.0.2 >R x64 3.2.3 （上のバージョンだと古すぎるので移行） これに従う． https://github.com/stan-dev/rstan/wiki/RStan-Getting-Started （１）基本の整備 バージョンにあったRtoolsのインストー…

RのMCMCライブラリのRstanの導入 windows 8 64bit R x64 3.0.2 これに従う． https://github.com/stan-dev/rstan/wiki/RStan-Getting-Started （１）基本の整備 バージョンにあったRtoolsのインストール （２）gccのパスをRtoolsのものにする（３）インスト…

RのMCMCライブラリのRstanの導入 windows 8 64bit R x64 3.0.2 これに従う． https://github.com/stan-dev/rstan/wiki/RStan-Getting-Started （１）基本の整備 バージョンにあったRtoolsのインストール （２）gccのパスをRtoolsのものにする（３）インスト…