【R言語】統計解析フリーソフトR 第6章【GNU R】 [無断転載禁止]©2ch.net

1132人目の素数さん2017/08/03(木) 19:23:12.67ID:Hq1blL0O
R は統計計算とグラフィックスのための言語・環境です。
統計計算で重宝するデータ型や、複数要素を処理する演算や関数、
解析結果を表示するグラフィックなど、多彩な機能を提供します。

●関連サイト
The R Project
http://www.r-project.org/
RjpWiki
http://www.okada.jp.org/RWiki/
リンク集
http://www.okada.jp.org/RWiki/?%A5%EA%A5%F3%A5%AF%BD%B8
※前スレ
【R言語】統計解析フリーソフトR 第5章【GNU R】
http://rio2016.2ch.net/test/read.cgi/math/1380168442/

525132人目の素数さん2018/09/11(火) 09:08:32.20ID:QUqp/jpE
数値を非零かどうか論理値に変換して処理するかは
関数によるな。

if やwhileは変換処理しているけど
whichだとエラーになるな。

> which(!0)
[1] 1
> which(0)
Error in which(0) : argument to 'which' is not logical
> which(2)
Error in which(2) : argument to 'which' is not logical
> which(!3)
integer(0)
>

526132人目の素数さん2018/09/11(火) 09:12:56.04ID:QUqp/jpE
!ってベクトル対応しているな
> which(!(-1:1))
[1] 2
> !(-2:2)
[1] FALSE FALSE TRUE FALSE FALSE

527132人目の素数さん2018/09/15(土) 08:05:49.52ID:Vl7XZ52q
!をnotではなくis.falseとかis.zero変換すれば( ・∀・)イイ!!

528132人目の素数さん2018/09/15(土) 16:42:00.50ID:h0gUCZ3r
>>527
is.zero だと?
だったら最初から == 0 てよいじゃないか。

529132人目の素数さん2018/09/15(土) 16:45:38.59ID:h0gUCZ3r
>>525
だいたい、if や while は関数じゃないし。
while (1) などと書くのはCに毒されているんじゃねーの?
Rなら while (TRUE) のほうが良いし、repeat というのもある。

530132人目の素数さん2018/09/15(土) 16:45:41.31ID:Vl7XZ52q
>>528
スクリプトを読むときに脳内変換する癖をつけておけば
可読性が悪いと思わずにすむ。

531132人目の素数さん2018/09/15(土) 16:46:13.42ID:Vl7XZ52q
>>518
>521みた?
if(length(x)!=0)になってないよ。

532132人目の素数さん2018/09/15(土) 17:19:35.45ID:h0gUCZ3r
>>531
Rのライブラリ?の書き方なんてあてにならないよ。
そもそも、rev も 509 の reverse も演算の回数を必要最低限にするなら、条件は length(x)) >= 2 または > 1 だ。
まあ if (length(x)) のほうが速かったのかも知れないが、常にそうなるとは限らないのだから、早すぎる最適化は諸悪の根源だ。

533132人目の素数さん2018/09/15(土) 19:30:34.55ID:Vl7XZ52q
可読性の次は速度かよ。
自分の流儀と違う { の使い方でも文句いいそう。

534132人目の素数さん2018/09/15(土) 19:44:27.43ID:VibLIqgl
0,1を1000万個作って
!
!=
>
==
での真偽判定を各々10回する時間を出してみた。


> x=rbinom(1e7,1,0.5)
> system.time(replicate(10,!x))
user system elapsed
1.27 0.33 1.59
> system.time(replicate(10,x!=0))
user system elapsed
2.39 0.36 2.75
> system.time(replicate(10,x>0))
user system elapsed
2.58 0.31 2.92
> system.time(replicate(10,x==1))
user system elapsed
2.60 0.36 2.99

535132人目の素数さん2018/09/15(土) 19:57:45.72ID:VibLIqgl
!x と x!=0 で10回やったが、

> f <- function(){
+ x=rbinom(1e7,1,0.5)
+ (a=system.time(!x)[3])
+ (b=system.time(x!=0)[3])
+ a<b
+ }
> (re=replicate(10,f()))
elapsed elapsed elapsed elapsed elapsed elapsed elapsed elapsed elapsed
TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
elapsed
TRUE

結果は再現された。

536132人目の素数さん2018/09/15(土) 20:06:37.73ID:VibLIqgl
lengthが絡むと

> x=1:1e8
> f1 = function(x) if(!length(x)) return(NULL) else x=x[-1]
> f2 = function(x) if(length(x)==0) return(NULL) else x=x[-1]
> system.time(f1(x))
user system elapsed
2.38 0.53 2.93
> system.time(f2(x))
user system elapsed
2.05 0.61 2.69

御指摘の通り、!の負け

537132人目の素数さん2018/10/04(木) 17:02:48.04ID:247Ted9r
>>444
>遅レスだが、Ubuntu 18.04LTSならデフォルトのIBus-mozcで日本語入力ができるそうだ。
>詳しくはUbuntu Weekly Recipeの第527回を読んでくだされ

Ubuntu 18.04.1 LTSにアップグレードしたけどできない
いままでスクレイピングにつかってrvestのログインコマンドもエラー出るようになったし最悪だ。。。。

538132人目の素数さん2018/10/04(木) 20:04:49.68ID:W6fWyA5T
アップグレードってことは文字入力がFcitxのままになってない?Fcitxならパッチ要だよ。
Voyager 18.04LTS(xubuntu 18.04LTS)ではパッチなしでRStudioに日本語入力できた。

539132人目の素数さん2018/10/26(金) 12:58:59.43ID:8LsI2NoU
昨日、PCを新しくして最新のR入れたんだけど、
作業ディレクトリの固定化できなくね?
前まではショートカットのプロパティで作業フォルダ指定してたんだけど、
新しいR、何度やっても作業ディレクトリがデフォルトのドキュメントに戻りやがる

どうすればいいか、誰がご教示を。

540132人目の素数さん2018/10/26(金) 13:04:50.23ID:I9GOHEF6
>>539
プロパティでコマンドラインオプションを消せばよい。
--cd-ほにゃらら
というのが付いてるでしょ。

541132人目の素数さん2018/10/26(金) 13:13:30.11ID:8LsI2NoU
>>540
あ、なるほど
ありがとうございます、多謝

542132人目の素数さん2018/10/26(金) 23:05:34.11ID:5yNEfvRx
xtsの取り扱い方について質問失礼します。
以下のコードで、dataには2018/10/20 00:00〜2018/10/25 00:00の1時間刻みのaとbの値のxts型オブジェクトが入りますが、(aとbは時系列データのイメージです)
このdataの毎日09:00の値を添え字を使って取り出す方法はありますか?
(data["2018-10-21 00:00"]とすると10/21の00:00が取り出せるような感じで)

library(xts)
time <- seq(as.POSIXct("2018-10-20 00:00"), as.POSIXct("2018-10-25 00:00"), by="1 hour")
a <- rnorm(length(time))
b <- rnorm(length(time))
data <- xts(x = cbind(a,b), order.by = time)

思いついたのはfor文を使って以下のようにするのかなと思ったのですがもう少しいい方法があるだろうなと思いまして…

data_result <- data[10] # まずxts型の入れ物を作る(2018/10/20 09:00)
for (i in 1:4) {
data_result <- rbind(temp,data[10+24*i]) # 2日目以降はこれにくっつける
}
data_result

543132人目の素数さん2018/10/26(金) 23:12:48.52ID:Xg1/ChR5
失礼、>>542の下のコードはこっちです

data_result <- data[10] # まずxts型の入れ物を作る(2018/10/20 09:00)
for (i in 1:4) {
data_result <- rbind(data_result, ,data[10+24*i]) # 2日目以降はこれにくっつける
}
data_result

544132人目の素数さん2018/10/26(金) 23:38:29.60ID:oEWqNGH8
>>543
idx <- seq(10, length(time), by=24)
result <- data[idx]
じゃダメなのか?

545132人目の素数さん2018/10/27(土) 01:10:34.48ID:NQVW6zPX
>>544
ああ、これがいいですね!お恥ずかしい…
ありがとうございます

546132人目の素数さん2018/11/01(木) 22:44:27.08ID:xCdOvDq8
巨大数を扱えるというふれこみのRmpfrって正確じゃないな。
50の階乗を計算させてみた

R with Rmpfr
> mpfr(factorial(50),1e5)
1 'mpfr' number of precision 100000 bits
30414093201713018969967457666435945132957882063457991132016803840

Haskell:
Prelude> product[1..50]
30414093201713378043612608166064768844377641568960512000000000000

Python:
import math
print(math.factorial(50))
30414093201713378043612608166064768844377641568960512000000000000

Wolfram:
https://www.wolframalpha.com/input/?i=50!
30414093201713378043612608166064768844377641568960512000000000000

547132人目の素数さん2018/11/02(金) 23:29:01.53ID:fMme23mF
こまけぇこたぁいいんだよ!!

548132人目の素数さん2018/11/03(土) 09:32:49.16ID:NiwdVATo
R3.5.xから、Windows環境だと日本語パスまわりでエラーがでてどうにもならない
enc2native()しても駄目だし
とりあえずcairo_png()だけでも…

549132人目の素数さん2018/11/03(土) 10:31:12.42ID:pqZePX+I
一度つかって日本語周りの処理でエラー出てきたんで未だに3.4.4浸かってる

550132人目の素数さん2018/11/04(日) 21:55:16.95ID:lTCeMsqQ
統計の質問はこのスレでいいんだろうか?

ある医院に1時間あたり平均5人の患者が来院し、その人数の分布はポアソン分布にしたがうとする。
1時間あたりの平均診療人数は6人で、一人あたりの診療時間は指数分布に従うとする。
診察までの平均の待ち時間は何時間か?

このシミュレーション解はこれであってる?

N=1e6
λ=5
μ=6
sum(rpois(N,λ)*rexp(N,μ))/N

> sum(rpois(N,λ)*rexp(N,μ))/N
[1] 0.8331036

60かけて、分にすると
[1] 49.86565

551132人目の素数さん2018/11/08(木) 12:35:54.73ID:wKTjJ6Fa
>>306
grep 使って書いてみた。
# mhs(c(1,0,1,1,0,1,1,1)) : return 3
mhs = function(x){ # maximum head sequence
y=paste(x,collapse='')
str="1"
if(!grepl(str,y)) return(0)
else{
while(grepl(str,y)){
str=paste0(str,"1")
}
return(nchar(str)-1)
}
}

system.time(mean(replicate(1e4,mhs(sample(0:1,100,rep=T))>=5)))


# N(=100)回コインをn(=5回)以上続けて表がでるか?TRUE or FALSE
seqn<-function(n=5,N=100,p=0.5){
rn=rbinom(N,1,p)
count=0
for(i in 1:N){
if(rn[i] & count<n){
count=count+1
}
else{
if(count==n) {return(TRUE)}
else{
count=0
}
}
}
return(count==n)
}

system.time(mean(replicate(1e4,seqn())))


結果は、
> system.time(mean(replicate(1e4,mhs(sample(0:1,100,rep=T))>=5)))
user system elapsed
4.56 0.01 4.61
> system.time(mean(replicate(1e4,seqn())))
user system elapsed
1.05 0.00 1.07

逆に4倍時間がかかるようになった

552132人目の素数さん2018/11/08(木) 14:05:00.34ID:1h1GfW+d
今さらだけどそれってrle使っちゃだめなの?

553132人目の素数さん2018/11/08(木) 15:13:38.19ID:UNhl34kg
coinの表裏を1と0で表して、その累積和を取ったベクトルを用意して
そのベクトルの5個前とのdiffの要素の中に5が一回でも現れることが、一回でも5回連続で表が出たことがあることの必要十分条件
grepはどうしても時間がかかるしif文もできれば使いたくない

searchseq <- function(n=5,N=100,p=0.5,trial=10000){
result <- 0 # 表が5回以上出た回数を数える入れ物
for (i in 1:trial){
coin <- rbinom(N,1,p)
coincumsum <- cumsum(coin)
coindiff <- diff(coincumsum,5)
#diff(coincumsum,5)の要素に一個でも5があれば表が5回以上出たことがあったということ
#anyでT/Fにして、sumで0/1にする
result <- result+sum(any(coindiff==5))
}
return(result/trial)
}

結果もよさそう
> searchseq()
[1] 0.80793

> system.time(mean(replicate(10000,mhs(sample(0:1,100,rep=T))>=5)))
ユーザ システム 経過
0.92 0.00 0.93
> system.time(mean(replicate(10000,seqn())))
ユーザ システム 経過
0.29 0.00 0.30
> system.time(searchseq())
ユーザ システム 経過
0.21 0.00 0.20

554132人目の素数さん2018/11/08(木) 22:10:55.78ID:I6htNxzg
あっ、コード中の5はnに置き換えてね

555132人目の素数さん2018/11/08(木) 23:32:48.68ID:/ZlhhGjJ
>>553
レスありがとうございました。

rleとも比べてみました。

> system.time(mean(replicate(1e4,seqn())))
user system elapsed
4.290 0.000 4.377
>system.time(mean(replicate(1e4,max(rle(rbinom(100,1,0.5))$len)>=5)))
user system elapsed
3.620 0.000 3.742
> system.time(mean(replicate(1e4,mhs(rbinom(100,1,0.5)>=5))))
user system elapsed
2.390 0.000 2.435
> system.time(searchseq())
user system elapsed
1.880 0.000 1.988

556132人目の素数さん2018/11/08(木) 23:35:10.38ID:/ZlhhGjJ
>>555
greplの逆転はsampleをrbinomに変えたためでしょう。

> system.time(replicate(1e5,sample(0:1,100,rep=T)))
user system elapsed
7.200 0.180 7.591
> system.time(replicate(1e5,rbinom(100,1,0.5)))
user system elapsed
5.980 0.230 6.319

557132人目の素数さん2018/11/09(金) 00:20:12.33ID:BZRFS9lT
不勉強ながらrle関数って初めて知ったけど使いやすそうだな

558132人目の素数さん2018/11/09(金) 00:36:46.75ID:Ui6BpICy
>>555
rleは0も数えているから間違っているんじゃ?

559132人目の素数さん2018/11/09(金) 02:01:36.95ID:Qla0VxTD
>>558
ご指摘ありがとうございました。
修正しました。
rle1 <- function(N=100,n=5,p=0.5){
r=rle(rbinom(N,1,p))
max(r$len[which(r$val==1)])
}
結果
> system.time(mean(replicate(1e4,max(rle1()>=5))))
user system elapsed
4.430 0.000 4.546
> system.time(mean(replicate(1e4,seqn())))
user system elapsed
4.490 0.010 4.609
> system.time(mean(replicate(1e4,mhs(rbinom(100,1,0.5))>=5)))
user system elapsed
7.440 0.000 7.656
> system.time(searchseq())
user system elapsed
1.950 0.000 2.066
>

560132人目の素数さん2018/11/09(金) 07:25:00.97ID:Qla0VxTD
無理矢理1行にして実行

system.time(mean(replicate(1e4,any(diff(cumsum(rbinom(100,1,0.5)),5)==5))))
user system elapsed
1.820 0.000 1.886
>
> system.time(mean(replicate(1e4,with(rle(rbinom(100,1,0.5)), max(lengths[wh
<e(1e4,with(rle(rbinom(100,1,0.5)), max(lengths[whi ch(values==1)])>=5))))
user system elapsed
4.370 0.010 4.478

561132人目の素数さん2018/11/09(金) 07:59:47.47ID:rijPDFSR
>>560
意味もなくforループ回してた上に毎回sum使って真偽値を数値に変換してたけど
replicate使って最後に一回だけmean取ると2.066→1.886で1割短くなるのね
他人のコード読むのは勉強になる

562132人目の素数さん2018/11/09(金) 09:06:56.99ID:Ui6BpICy
>>559
whichいらねーよ。

563132人目の素数さん2018/11/09(金) 09:44:07.41ID:5AnUTlVm
>>559
全部が0のとき、エラーになるので修正

rle01 <- function(x){ # c(0,1,1,1,0,0) => return 3
if(sum(x)==0) return(0) #c(0,0,0,0,0,0) => return 0
else{
r=rle(x) # Run Length Encoding
max(r$lengths[which(r$values==1)]) # max length of value 1
}
}

動作確認

> rle01(x<-rbinom(100,1,0.5)) ; x
[1] 8
[1] 1 1 0 1 0 1 0 0 0 1 1 1 0 1 0 0 0 0 1 0 0 0 1 1 1 1 0 0 1 1 1 1 1 0 0
[36] 1 1 0 1 1 1 0 1 1 0 0 1 1 0 1 1 1 0 1 1 0 0 1 1 0 1 0 1 1 0 1 0 1 0 0
[71] 1 1 1 1 1 1 0 1 1 1 0 1 1 1 1 1 1 1 1 0 0 1 0 0 1 1 1 0 1 1
> rle01(x<-rbinom(100,1,0)) ; x
[1] 0
[1] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
[36] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
[71] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0

564132人目の素数さん2018/11/09(金) 10:10:01.39ID:Ui6BpICy
>>563
センスないなあ。
rle01 <- function(x) {
r <- rle(x)
one <- r$values == 1
if (any(one)) max(r$length[one]) else 0
}

565132人目の素数さん2018/11/09(金) 10:40:31.84ID:5AnUTlVm
>>564
whichは余分でした。

sumの方がrleより高速だと思ったらから、すべて0の場合はrleを呼ばないことにしただけ。
0連続でやるとこれだけ差がつく。

rle01 <- function(x){ # c(0,1,1,1,0,0) => return 3
if(sum(x)==0) return(0) #c(0,0,0,0,0,0) => return 0
else{
r=rle(x) # Run Length Encoding
max(r$lengths[which(r$values==1])] # max length of value 1
}
}
rle012 <- function(x) {
r <- rle(x)
one <- r$values == 1
if (any(one)) max(r$length[one]) else 0
}

> x=rep(0,1e8)
> system.time(rle01(x))
user system elapsed
0.3 0.0 0.3
> system.time(rle012(x))
user system elapsed
7.36 4.52 13.25
>

566132人目の素数さん2018/11/09(金) 10:52:23.39ID:Ui6BpICy
>>565
そんな特別な場合のことに対して高速化するのは愚の骨頂。
1が一つでもある場合にsumを呼ぶのは余計じゃないのか?
余計なwhichがあるくらいだから、まずは素直にやりたいこと、やるべきことを正しく書くようにしたら?

567132人目の素数さん2018/11/09(金) 11:01:40.74ID:5AnUTlVm
1000本に1本あたる宝クジを100本買って続けて2本あたる確率のシミュレーション解の算出時間比較。
確率の理論値は9.8897353347449091e-05

> system.time(mean(replicate(1e4,rle01(rbinom(N,1,p))>=n)))
user system elapsed
0.61 0.00 0.64

> system.time(mean(replicate(1e4,rle12(rbinom(N,1,p))>=n)))
user system elapsed
1.97 0.00 2.03

568132人目の素数さん2018/11/09(金) 15:22:38.24ID:5AnUTlVm
>>567
分数で表すと

890788167367/9007199254740992

9.889735334744909e-05

569132人目の素数さん2018/11/09(金) 19:43:44.34ID:0M0agNOf
JPXのデータ、ファイル形式csvを読み込もうとするとうまく行かないんですが
どんな引数をつければいいですか

570132人目の素数さん2018/11/10(土) 10:36:49.41ID:QJ6NJqU7
>>566
実はシミュレーションじゃなくて
漸化式からのプログラム解を分数表示するプログラムはpythonで作成済。
ここに置いた。
https://egg.2ch.net/test/read.cgi/hosp/1540905566/77

571132人目の素数さん2018/11/10(土) 11:20:41.80ID:QJ6NJqU7
コインを100回ふったときの表連続の最大数が5であったときの
このコインの表がでる確率の期待値、モード比、信頼区間を求めるのが次のネタ。

unirootで算出できたけど
シミュレーションはどうすればいいのかアイデアが浮かばない。
MCMCで解決できるかなぁ?

これをシミュレーションで検証したい。

$Rscript main.r
lower mean mode upper
0.2487456 0.4469764 0.4589692 0.6386493

http://tpcg.io/asKRE9

572132人目の素数さん2018/11/11(日) 20:11:21.81ID:ODPKEEGK
1億回のコイントスで何回連続して表がでる確率が高いかRでやってみた。

# maximal sequential head probability at 10^8 coin flip
> y
[1] 2.2204460492503131e-16 2.2204460492503131e-16
[3] 8.8817841970012523e-16 1.5543122344752192e-15
[5] 3.5527136788005009e-15 6.8833827526759706e-15
[7] 1.4210854715202004e-14 2.8199664825478976e-14
[9] 5.6843418860808015e-14 1.1346479311669100e-13
[11] 2.2737367544323206e-13 4.5452530628153909e-13
[13] 9.0949470177292824e-13 1.8187673589409314e-12
[15] 3.6379788070917130e-12 7.2757355695785009e-12
[17] 1.4551915228366852e-11 2.9103608412128779e-11
[19] 5.8207660913467407e-11 1.1641509978232989e-10
[21] 6.6493359939245877e-06 2.5720460687386204e-03
[23] 4.8202266324911647e-02 1.7456547460031679e-01
[25] 2.4936031630254019e-01 2.1428293501123058e-01
[27] 1.4106434838399229e-01 8.1018980443629832e-02
[29] 4.3428433624081136e-02 2.2484450838189007e-02

25回が続くのが4回に1回あることになる。
pythonで25回以上と25回ちょうどになるのを計算させてみた。

その結果、
Over 25
6977459029519597/9007199254740992
= 0.7746535667951356
Just 25
2246038053550679/9007199254740992
= 0.24936031612362342

高速化を狙ってCに移植したら
100万回で暴走。
https://egg.5ch.net/test/read.cgi/hosp/1540905566/132

573132人目の素数さん2018/11/12(月) 21:04:57.19ID:boi8bvdM
"マラソン大会の選手に1から順番に番号の書かれたゼッケンをつける。
用意されたゼッケンN(=100)枚以下の参加であった。
無作為に抽出したM(=5)人のゼッケン番号の最大値はMmax(=60)であった。
参加人数推定値の期待値とその95%信頼区間を求めよ"

decken <- function(M, Mmax, N, conf.level=0.95){
if(Mmax < M) return(0)
n=Mmax:N
pmf=choose(Mmax-1,M-1)/choose(n,M)
pdf=pmf/sum (pmf)
mean=sum(n*pdf)
upr=n[which(cumsum(pdf)>conf.level)[1]]
lwr=Mmax
c(lower=lwr,mean=mean,upper=upr)
}
decken(M=5,Mmax=60,N=100)

> decken(M=5,Mmax=60,N=100)
lower mean upper
60.0000 71.4885 93.0000

これをシミュレーションで確認したい。

# simulation
M=5 ; Mmax=60 ; N=100
sub <- function(M,Mmax,N){
n=sample(Mmax:N,1) # n : 参加人数n
m.max=max(sample(n,M)) # m.max : n人からM人選んだ最大番号
if(m.max==Mmax) return(n)
}
sim <- function(){
n=sub(M,Mmax,N)
while(is.null(n)){
n=sub(M,Mmax,N) # 最大番号が一致するまで繰り返す
}
return(n)
}
runner=replicate(1e4,sim())
summary(runner) ; hist(runner,freq=F,col="lightblue")
quantile(runner,prob=c(0,0.95))
cat(names(table(runner)[which.max(table(runner))]))

> summary(runner) ; hist(runner,freq=F,col="lightblue")
Min. 1st Qu. Median Mean 3rd Qu. Max.
60.00 63.00 68.00 71.43 77.00 100.00
> quantile(runner,prob=c(0,0.95))
0% 95%
60 93
> cat(names(table(runner)[which.max(table(runner))])) # 最頻値
60

結果は確認できたけど、もっと高速なシミュレーションアルゴリズムはあるだろうか?

574132人目の素数さん2018/11/16(金) 13:44:59.80ID:U19cHKqd
重複があるか否かを返す、anyDuplicatedという関数を知ったので総当たり比較と早いかどうか比べてみた。
覆面算を □/□ * □/□ = □□/□を解くの使ってみた。
a/b * c/d == ef/g (c>a)として

F = function(fun){
n=1:9
ans=NULL
for(a in n){
for(b in n){
for(c in a:9){
for(d in n){
for(e in n){
for(f in n){
for(g in n){
if(fun(a,b,c,d,e,f,g)){
ans=rbind(ans,c(a,b,c,d,e,f,g))}}}}}}}}

return(ans)
}
で虱潰しに判定

F1=function(a,b,c,d,e,f,g){ #全部の組み合わせが等しくないのを確認
(a/b)*(c/d)==(10*e+f)/g &
a!=b & a!=c & a!=d & a!=e & a!=f & a!=g &
b!=c & b!=d & b!=e & b!=f & b!=g & c!=d &
c!=e & c!=f & c!=g & d!=e & d!=f & d!=g & e!=f & e!=g & f!=g
}

F2=function(a,b,c,d,e,f,g){ # anyDuplicatedで重複なしを判定
(a/b)*(c/d)==(10*e+f)/g & !anyDuplicated(c(a,b,c,d,e,f,g))
}

> system.time(F(F1))
user system elapsed
52.56 0.25 53.38
> system.time(F(F2))
user system elapsed
113.78 0.11 115.81

anyDuplicatedでコードは短くなるが速さが犠牲になった。

575132人目の素数さん2018/11/16(金) 16:13:49.03ID:zcEc4+wx
カルマンフィールタ

新着レスの表示
レスを投稿する