17/10/11 22:25:55.39 HJHHXjdz.net
>>155
6行4列の左上の6つ1,2,3,7,8,9の平均5だと思う。
158:１３２人目の素数さん
17/10/11 22:59:20.81 HJHHXjdz.net
>>154
効率的かどうかはわからんが、やってみた
(d=array(1:48, dim=c(6, 4, 2)))
a=3
b=2
l
159:=6 m=4 n=2 re=NULL for(k in 1:n){ for(j in 1:(m/b)){ for(i in 1:(l/a)){ re=append(re,mean(d[1:3+3*(i-1),1:2+2*(j-1),k])) } } } array(re,dim=c(m/b,l/a,n))
160:１３２人目の素数さん
17/10/11 23:00:23.97 UQm0Ph1Q.net
>>155
>>154です。
説明が悪くてすみません。
>>156の通りです。
161:１３２人目の素数さん
17/10/11 23:02:02.80 UQm0Ph1Q.net
>>157
ありがとうございます。
手元にデータがないので明日試してみます。
162:１３２人目の素数さん
17/10/11 23:02:32.81 bFq2a2U4.net
なるほど、じゃあその5,8,17,20はこういうことか。
> a <- array(1:48, dim=c(6, 4, 2))
> b <- expand.grid(1:(dim(a)[1] %/% 3), 1:(dim(a)[2] %/% 2))
> f <- function(i, j, k=1){mean(a[(3*i-2):(3*i),(2*j-1):(2*j),k])}
> mapply(f, b[,1], b[,2])
[1] 5 8 17 20
163:１３２人目の素数さん
17/10/11 23:04:06.69 bFq2a2U4.net
あぁ、すでに >>157 が�
164:嘯曹｢ていたorz
165:１３２人目の素数さん
17/10/11 23:06:59.19 HJHHXjdz.net
> d=array(1:48, dim=c(6, 4, 2))
> a=3
> b=2
> l=6
> m=4
> n=2
> re=NULL
> for(k in 1:n){
+ for(j in 1:(m/b)){
+ for(i in 1:(l/a)){
+ re=append(re,mean(d[1:a+a*(i-1),1:b+b*(j-1),k]))
+ }
+ }
+ }
> array(re,dim=c(m/b,l/a,n))
, , 1
[,1] [,2]
[1,] 5 17
[2,] 8 20
, , 2
[,1] [,2]
[1,] 29 41
[2,] 32 44
データは10000x10000x400なら
l=10000
m=10000
n=400
でよくね？
166:１３２人目の素数さん
17/10/12 07:10:38.11 2sFs7Gsb.net
最初の1行で
> d=array(1:(10000*10000*400))
Error: cannot allocate vector of size 298.0 Gb
エラーがでたｗ
167:１３２人目の素数さん
17/10/12 07:12:09.84 2sFs7Gsb.net
正しいエラーｗはこっち
> d=array(1:(10000*10000*400),dim=c(6, 4, 2))
Error: cannot allocate vector of size 298.0 Gb
168:１３２人目の素数さん
17/10/12 07:21:02.01 2sFs7Gsb.net
そもそも、10000x10000x400だと10000行は3行にならないな。最終行を無視して9999行にするのかな？
169:１３２人目の素数さん
17/10/12 14:48:46.91 sG1WSodN.net
>160を参考に>157のforだらけのスパゲティー・ソースを（ちょっと一般化して）改変しました。
（変数名の重複を避けるため>160の変数名は変更してあります。）
> a=3
> b=2
> c=1
>
> l=6
> m=4
> n=2
>
> A <- array(1:48, dim=c(l,m,n))
> B <- expand.grid(1:(l %/% a), 1:(m %/% b),1:(n %/% c))
> f <- function(i, j, k){mean(A[(a*(i-1)+1):(a*i),(b*(j-1)+1):(b*j),(c*(k-1)+1):(c*k)])}
> re=mapply(f, B[,1], B[,2],B[,3])
> array(re,dim=c(m%/%b,l%/%a,n%/%c))
, , 1
[,1] [,2]
[1,] 5 17
[2,] 8 20
, , 2
[,1] [,2]
[1,] 29 41
[2,] 32 44
>
上級者のコードを読むのは勉強になります。（深謝）
170:１３２人目の素数さん
17/10/12 16:33:51.03 f/ETd6de.net
先輩方，>>154です。
巨細にわたってご教示・ご議論いただきありがとうございました。
>>166のコードを参考にして，無事処理することができました。
実行速度も申し分なく，効率的に作業が進められそうです。
今回は私には目から鱗な解法で，勉強になりました。
質問して良かったと思っています。。本当にありがとうございました。
10,000行の扱いですが，3行ごとに処理して1行無視してもよいし，
4から7行程度で処理してもよく，フレキシブルです。
今回は4行にして余りなしにしました。
171:１３２人目の素数さん
17/10/12 18:43:54.51 2sFs7Gsb.net
処理速度を比べてみた。
fs<-function(){ #スパゲッティ・ソース
d=array(1:(l*m*n),dim=c(l,m,n))
re=NULL
for(k in 1:(n/c)){
for(j in 1:(m/b)){
for(i in 1:(l/a)){
re=append(re,mean(d[1:a+a*(i-1),1:b+b*(j-1),1:c+c*(k-1)]))
}}}
array(re,dim=c(m/b,l/a,n))
}
fe <- function(){ # Expertソース
A <- array(1:(l*m*n), dim=c(l,m,n))
B <- expand.grid(1:(l %/% a), 1:(m %/% b),1:(n %/% c))
f <- function(i, j, k){mean(A[(a*(i-1)+1):(a*i),(b*(j-1)+1):(b*j),(c*(k-1)+1):(c*k)])}
re=mapply(f, B[,1], B[,2],B[,3])
array(re,dim=c(m%/%b,l%/%a,n%/%c))
}
a=4
b=2
c=1
l=400
m=200
n=2
> system.time(fs())
user system elapsed
2.17 0.00 2.18
> system.time(fe())
user system elapsed
1.14 0.00 1.14
>160氏のスキルに改めて感服。
172:１３２人目の素数さん
17/10/14 23:02:55.66 WIgy7a2M.net
# 既約剰余類群 (Z/nZ)*
kjrg <- function(n){
nn=0:(n-1)
f=function(x,y) (x*y)%%n
names(nn)=paste0('n',nn)
z=outer(nn,nn,f)
idx=which(gmp::gcd(1:(n-1),n)==1)+1
return(z[idx,idx])
}
kjrg(10)
kjrg(15)
kjrg(24)
> kjrg(24)
n1 n5 n7 n11 n13 n17 n19 n23
n1 1 5 7 11 13 17 19 23
n5 5 1 11 7 17 13 23 19
n7 7 11 1 5 19 23 13 17
n11 11 7 5 1 23 19 17 13
n13 13 17 19 23 1 5 7 11
n17 17 13 23 19 5 1 11 7
n19 19 23 13 17 7 11 1 5
n23 23 19 17 13 11 7 5 1
173:１３２人目の素数さん
17/10/15 13:21:32.63 GlVyhxtn.net
# 既約剰余類群での加減乗除
n=24
n_star=which(gmp::gcd(1:(n-1),n)==1) ; n_star
tasu <- function(x,y) (x+y)%%n
hiku <- function(x,y) (x-y)%%n # row - col
kake <- function(x,y) (x*y)%%n
g=function(x) n_star[which(x==
174:1)] .M=outer(n_star,n_star,kake) G=apply(.M,2,g) gyaku <- function(x) n_star[which(G==(x%%n))] waru <- function(x,y) (x*gyaku(y))%%n # row ÷ col xx=yy=c(0,n_star) names(xx)=paste0('x',c(0,n_star)) names(yy)=paste0('y',c(0,n_star)) outer(xx,yy,tasu) # x + y outer(xx,yy,hiku) # x - y outer(xx,yy,kake) # x * y X=Y=n_star #outer(X,Y,waru) # WRONG!! a=expand.grid(X,Y) b=matrix(mapply(waru,a[,1],a[,2]),length(X)) rownames(b)=paste0('x',n_star) colnames(b)=paste0('y',n_star) b # x / y
175:１３２人目の素数さん
17/10/15 14:03:02.02 GlVyhxtn.net
名前も一致してidenticalなんだな。
> x=1:10
> y=1:10
> identical(x,y)
[1] TRUE
> names(x)=LETTERS[1:10]
> identical(x,y)
[1] FALSE
> names(y)=LETTERS[1:10]
> identical(x,y)
[1] TRUE
> names(y)=letters[1:10]
> identical(x,y)
[1] FALSE
allの方は内容の全体が一致していれば順番は入れ替わってもTRUEなんだな。
> x=1:10
> y=1:10
> all(x,y)
[1] TRUE
> names(x)=LETTERS[1:10]
> all(x,y)
[1] TRUE
> names(y)=letters[1:10]
> all(x,y)
[1] TRUE
> all(x,sample(y))
[1] TRUE
176:１３２人目の素数さん
17/10/16 01:11:00.56 NsB0m+v3.net
beki <- function(x,n,p){ # x^n%%p
if(n==0) return(1)
if(n==1) return(x%%p)
re=numeric(n)
re[1]=x
for(i in 1:(n-1)) re[i+1] = (x*re[i])%%p
return(re[n])
}
p=1009 # フェルマーの小定理
xx=0:(p-1)
n=p-1
y=sapply(xx,function(x)beki(x,n,p))
summary(y)
177:１３２人目の素数さん
17/10/16 07:46:39.05 NsB0m+v3.net
90歳女性のお看取り
>120のコードから女性の平均余命をグラフにしてみた
Age=85:105
(Life_Expectancy=sapply(Age,function(x) LE(f,x)))
plot(Age,Life_Expectancy,pch=19,col=sample(colors()))
abline(h=5,col='gray',lty=2)
URLﾘﾝｸ(i.imgur.com)
平均余命が5年を切るのは92歳なので、老衰にせずに
主治医指定の病名で診断書を作成。
カルテ記載のしっかりした患者の診療は負担なくこなせて（・∀・）ｲｲ!!
178:１３２人目の素数さん
17/10/16 22:16:44.56 ACageh59.net
outerって使えない場面があるよね？
剰余系での加減乗除での実際
> n=5 # prime number
> nn=1:(n-1)
> tasu <- function(x,y) (x+y)%%n
> hiku <- function(x,y) (x-y)%%n # row - col
> kake <- function(x,y) (x*y)%%n
> g=function(x) nn[which(x==1)]
> .M=outer(nn,nn,kake)
> G=apply(.M,2,g)
> gyaku <- function(x) nn[which(G==(x%%n))]
> waru <- function(x,y) (x*gyaku(y))%%n # row / col
> waru(3,2)
[1] 4
> xx=yy=c(0,nn)
> names(xx)=paste0('x',c(0,nn))
> names(yy)=paste0('y',c(0,nn))
> outer(xx,yy,tasu) # x + y
y0 y1 y2 y3 y4
x0 0 1 2 3 4
x1 1 2 3 4 0
x2 2 3 4 0 1
x3 3 4 0 1 2
x4 4 0 1 2 3
> outer(xx,yy,hiku) # x - y
y0 y1 y2 y3 y4
x0 0 4 3 2 1
x1 1 0 4 3 2
x2 2 1 0 4 3
x3 3 2 1 0 4
x4 4 3 2 1 0
> outer(xx,yy,kake) # x * y
y0 y1 y2 y3 y4
x0 0 0 0 0 0
x1 0 1 2 3 4
x2 0 2 4 1 3
x3 0 3 1 4 2
x4 0 4 3 2 1
> X=Y=nn
> outer(X,Y,waru) # WRONG!!
[,1] [,2] [,3] [,4]
[1,] 1 1 1 1
[2,] NA NA NA NA
[3,] NA NA NA NA
[4,] NA NA NA NA
> a=expand.grid(X,Y)
> b=matrix(mapply(waru,a[,1],a[,2]),length(X))
> rownames(b)=paste0('x',nn)
> colnames(b)=paste0('y',nn)
> b # x / y
y1 y2 y3 y4
x1 1 3 2 4
x2 2 1 4 3
x3 3 4 1 2
x4 4 2 3 1
179:１３２人目の素数さん
17/10/16 23:28:56.68 FS7yMelZ.net
>>174
そりゃあ、outer に渡す関数がまずいんだろ。
outer に渡す関数は、完全にベクトル対応になっていなければならないんじゃないの？
でなければ、わざわざ outer 使うこともなかろう。
180:１３２人目の素数さん
17/10/17 06:58:05.54 ps5i5Wv8.net
>>175
レスありがとうございます。
自作関数から自作関数を呼び出しているのが原因かと思っていましたが、
ベクトル対応していないとはこういうことでしょうか？
加減乗は引数のどちらをベクトルにしても動作
> kake(1:4,1)
[1] 1 2 3 4
> kake(1,1:4)
[1] 1 2 3 4
> tasu(1:4,1)
[1] 2 3 4 0
> tasu(1,1:4)
[1] 2 3 4 0
> hiku(1:4,1)
[1] 0 1 2 3
> hiku(1,1:4)
[1] 0 4 3 2
除算は第２引数をベクトルにすると誤動作
> waru(1:4,1)
[1] 1 2 3 4
> waru(4,1:4)
[1] 4 1
181:１３２人目の素数さん
2017/10/1
182:7(火) 22:43:30.39 ID:yEir4ohi.net
183:名無しさん＠そうだ選挙に行こう！ Go to vote!
17/10/22 07:51:52.36 /ES4CbkI.net
>>177
ありがとうございます。
outerのソースをみたら、こんな記述があったので、X,Yを同時にベクトルで与えても動作する関数でないとダメみたいでした。
FUN <- match.fun(FUN)
Y <- rep(Y, rep.int(length(X), length(Y)))
if (length(X))
X <- rep(X, times = ceiling(length(Y)/length(X)))
FUN(X, Y, ...)
184:名無しさん＠そうだ選挙に行こう！ Go to vote!
17/10/22 08:28:19.35 /ES4CbkI.net
# greatest common divisor
gcd <- function(x,y) ifelse(x|y,ifelse(x*y,max(intersect((1:abs(x))[abs(x)%%(1:abs(x))==0],(1:abs(y))[abs(y)%%(1:abs(y))==0])),max(abs(x),abs(y))),NA)
非正整数対応すると長くなるなぁ
> gcd(24,15)
[1] 3
> gcd(16,15)
[1] 1
> gcd(5,0)
[1] 5
> gcd(-24,15)
[1] 3
> gcd(-24,-15)
[1] 3
> gcd(0,0)
[1] NA
185:名無しさん＠そうだ選挙に行こう！ Go to vote!
17/10/22 09:17:13.05 /ES4CbkI.net
ユークリッドの互除法をつかうとこんな感じ
> GCD <- function(x,y){
+ if(round(x)!=x | round(y)!=y) stop('Not integer')
+ if(!x&!y) return(NA)
+ a=max(c(abs(x),abs(y)))
+ b=min(c(abs(x),abs(y)))
+ if(!a*b) return(a)
+ r=integer()
+ r[1]=a
+ r[2]=b
+ i=1
+ r[i+2]=r[i]%%r[i+1]
+ while(r[i+2]!=0 & r[i+2]!=1){
+ i=i+1
+ r[i+2]=r[i]%%r[i+1]
+ }
+ return(ifelse(r[i+2]==0,r[i+1],1))
+ }
> GCD(0,0)
[1] NA
> GCD(24,15)
[1] 3
> GCD(16,15)
[1] 1
> GCD(5,0)
[1] 5
> GCD(24,-15)
[1] 3
> GCD(-24,-15)
[1] 3
> GCD(2.3,4)
Error in GCD(2.3, 4) : Not integer
186:名無しさん＠そうだ選挙に行こう！ Go to vote!
17/10/22 10:07:52.18 tFigTWxc.net
文系研究者のワイ、困惑
ディープな世界やでほんま
187:１３２人目の素数さん
17/10/25 20:38:52.07 A26rYIkF.net
#.n発r中の狙撃手が.N発狙撃するときの命中数を返す
Go13 <- function(.N, .n, r, k=10^3){ # k:シミュレーション回数
f <-function(S,N=.N,n=.n){
y=c(rep(1,S),rep(0,N-S))
sum(sample(y,n))
}
xx=r:.N
SS=NULL
for(i in 1:k){
x=sapply(xx,f)
SS=c(SS,which(x==r)-1+r)
}
print(summary(SS))
invisible(SS)
}
188:￥
17/10/26 04:11:35.13 jT09z118.net
￥
189:￥
17/10/26 04:11:58.98 jT09z118.net
￥
190:￥
17/10/26 04:12:18.59 jT09z118.net
￥
191:￥
17/10/26 04:12:35.94 jT09z118.net
￥
192:￥
17/10/26 04:12:53.44 jT09z118.net
￥
193:￥
17/10/26 04:13:10.90 jT09z118.net
￥
194:￥
17/10/26 04:13:28.27 jT09z118.net
￥
195:￥
17/10/26 04:13:46.40 jT09z118.net
￥
196:￥
17/10/26 04:14:03.59 jT09z118.net
￥
197:￥
17/10/26 04:14:22.57 jT09z118.net
￥
198:１３２人目の素数さん
17/10/26 18:17:25.90 tHqsRRzo.net
こういう問題を考えるのは楽しいね。
ゴルゴ13は100発100中
ゴルゴ14は10発10中
ゴルゴ15は1発1中
とする。
各々10000発撃ったとき各ゴルゴの命中数の期待値はいくらか？
ドツボ13は100発0中
ドツボ14は10発0中
ドツボ15は1発0中
とする。
各々10000発撃ったときドツボの命中数の期待値はいくらか？
199:１３２人目の素数さん
17/11/01 12:40:01.83 UKV+bz4p.net
ベクトルは再利用されるんだな
> c(1,2,1,2)==c(1,2)
[1] TRUE TRUE TRUE TRUE
> all(c(1,2,1,2)==c(1,2))
[1] TRUE
> equal <- function(x,y) sum((x-y)^2)==0 & length(x)==length(y)
> equal(c(1,2,1,2),c(1,2))
[1] FALSE
200:１３２人目の素数さん
17/11/03 08:27:39.47 7YwYHhPi.net
> x= list(length=5)
> x
$length
[1] 5
> y=list()
> length(y)=4
> y
[[1]]
NULL
[[2]]
NULL
[[3]]
NULL
[[4]]
NULL
> z=list(5)
> z
[[1]]
[1] 5
> w = vector(mode="list" , length= 5)
> w
[[1]]
NULL
[[2]]
NULL
[[3]]
NULL
[[4]]
NULL
[[5]]
NULL
201:１３２人目の素数さん
17/11/09 08:41:04.22 SxG1YX2j.net
Brunner-Munzel検定って比較する両群に重なりがないとｐ値がＮＡになるんだが、
これってバ
202:グなのか？ > x=1:5 > y=6:10 > lawstat::brunner.munzel.test(x,y) Brunner-Munzel Test data: x and y Brunner-Munzel Test Statistic = Inf, df = NaN, p-value = NA 95 percent confidence interval: NaN NaN sample estimates: P(X<Y)+.5*P(X=Y) 1 > x=1:6 > y=6:11 > lawstat::brunner.munzel.test(x,y) Brunner-Munzel Test data: x and y Brunner-Munzel Test Statistic = 24.749, df = 10, p-value = 2.651e-10 95 percent confidence interval: 0.9423463 1.0298759 sample estimates: P(X<Y)+.5*P(X=Y) 0.9861111
203:１３２人目の素数さん
17/11/14 22:02:36.50 Jlt4ZA5T.net
突然ですが質問があります。最近Rソフトを利用し始めた者のです。単刀直入に、
ボロノイ図を制作した場合に分割された面積の平均を出力するには、どうすればよいのでしょうか？
204:１３２人目の素数さん
17/11/14 22:06:32.42 Jlt4ZA5T.net
現状としては住所をグーグルマップの緯度経度に変換し、分割まではできています。
宜しければどなたかお力添えをいただけないでしょうか？
205:￥
17/11/15 00:58:17.28 WZuPK5Ir.net
￥
206:￥
17/11/15 00:58:38.82 WZuPK5Ir.net
￥
207:￥
17/11/15 00:58:56.67 WZuPK5Ir.net
￥
208:￥
17/11/15 00:59:14.81 WZuPK5Ir.net
￥
209:￥
17/11/15 00:59:32.78 WZuPK5Ir.net
￥
210:￥
17/11/15 00:59:50.57 WZuPK5Ir.net
￥
211:￥
17/11/15 01:00:09.89 WZuPK5Ir.net
￥
212:￥
17/11/15 01:00:29.29 WZuPK5Ir.net
￥
213:￥
17/11/15 01:00:47.81 WZuPK5Ir.net
￥
214:￥
17/11/15 01:01:09.56 WZuPK5Ir.net
￥
215:１３２人目の素数さん
17/11/15 17:11:14.46 eZCzfwkb.net
>>198
SpatialPolygonsクラスまでできるといると仮定して、
各ポリゴンの中にareaスロットがあるから、
それを抽出して、その平均を取れば良い。
ボロノイ図がSpatialPolygonsクラスになっていないなら、
まずは変換。
216:￥
17/11/15 20:30:49.23 WZuPK5Ir.net
￥
217:￥
17/11/15 20:31:06.31 WZuPK5Ir.net
￥
218:￥
17/11/15 20:31:21.95 WZuPK5Ir.net
￥
219:￥
17/11/15 20:31:40.54 WZuPK5Ir.net
￥
220:￥
17/11/15 20:31:58.15 WZuPK5Ir.net
￥
221:￥
17/11/15 20:32:16.33 WZuPK5Ir.net
￥
222:￥
17/11/15 20:32:32.65 WZuPK5Ir.net
￥
223:￥
17/11/15 20:32:48.32 WZuPK5Ir.net
￥
224:￥
17/11/15 20:33:05.48 WZuPK5Ir.net
￥
225:￥
17/11/15 20:33:21.41 WZuPK5Ir.net
￥
226:１３２人目の素数さん
17/11/15 21:12:29.48 9zgzmHZc.net
ボロノイ図の各領域の面積の平均値って、単に
描画領域の総面積÷点数 じゃないの？
何か暗黙の前提があるの？
227:￥
17/11/16 06:01:20.12 6ldUKvsQ.net
￥
228:￥
17/11/16 06:01:37.22 6ldUKvsQ.net
￥
229:￥
17/11/16 06:01:56.08 6ldUKvsQ.net
￥
230:￥
17/11/16 06:02:14.27 6ldUKvsQ.net
￥
231:￥
17/11/16 06:02:31.67 6ldUKvsQ.net
￥
232:￥
17/11/16 06:02:48.95 6ldUKvsQ.net
￥
233:￥
17/11/16 06:03:06.96 6ldUKvsQ.net
￥
234:￥
17/11/16 06:03:24.98 6ldUKvsQ.net
￥
235:￥
17/11/16 06:03:42.58 6ldUKvsQ.net
￥
236:￥
17/11/16 06:04:00.50 6ldUKvsQ.net
￥
237:１３２人目の素数さん
17/11/17 20:18:17.65 5lHMDqma.net
Stan って 離散量をparameter 扱いできないので
m ? categorical(θ) としても ｍ の分布を出せない。
JAGSだと
m ～ dcat(θ）
とすると ｍの分布がだせる。
この理解でいいかな？
238:￥
17/11/17 20:41:11.51 Y8c01QBt.net
￥
239:￥
17/11/17 20:41:27.83 Y8c01QBt.net
￥
240:￥
17/11/17 20:41:43.62 Y8c01QBt.net
￥
241:￥
17/11/17 20:41:58.17 Y8c01QBt.net
￥
242:￥
17/11/17 20:42:13.43 Y8c01QBt.net
￥
243:￥
17/11/17 20:42:28.14 Y8c01QBt.net
￥
244:￥
17/11/17 20:42:43.14 Y8c01QBt.net
￥
245:￥
17/11/17 20:42:59.36 Y8c01QBt.net
￥
246:￥
17/11/17 20:43:14.19 Y8c01QBt.net
￥
247:￥
17/11/17 20:43:30.82 Y8c01QBt.net
￥
248:１３２人目の素数さん
17/11/20 17:47:00.70 ptVXRSQM.net
以前、ボロノイ図について質問させていただいた初心者です。
質問後色々と試行錯誤したのですがdir.areaを出せず悩んでおります。
library(deldir)
> library(ggmap)
> df <- read.csv("C:/Users/***/Desktop/***/***.csv",header=T)
> loc <- c(min(df$lon),min(df$lat),max(df$lon),max(df$lat))
> GM <-ggmap(get_map(location=loc,zoom=11,source="google"))
> GP <- geom_point(data=df,aes(x=lon,y=lat,colour=factor(type)),size=3)
> dd <- deldir(df$lon,df$lat)
>GS <-geom_segment(data=dd$dirsgs,aes(x=x1,y=y1,xend=x2,yend=y2),size=0.5,linetype=1)
> GM+GP+GS
と入力しボロノイ図は出せているのですが、
249:１３２人目の素数さん
17/11/20 17:52:19.12 ptVXRSQM.net
summaryと入力しても
function (object, ...)
UseMethod("summary")
<bytecode: 0x06b5eef0>
<environment: namespace:base>
と表示されます。どの様に入力または再入力すればdir.areaが表示させられるでしょうか？
私自身初めてで、全く分からない状態でRを使用しているので、宜しければ分かり易く教えて頂けると幸いです。
250:￥
17/11/20 21:05:10.85 yjl62pX+.net
￥
251:￥
17/11/20 21:05:28.54 yjl62pX+.net
￥
252:￥
17/11/20 21:05:44.92 yjl62pX+.net
￥
253:￥
17/11/20 21:06:01.26 yjl62pX+.net
￥
254:￥
17/11/20 21:06:18.27 yjl62pX+.net
￥
255:￥
17/11/20 21:06:35.26 yjl62pX+.net
￥
256:￥
17/11/20 21:06:52.55 yjl62pX+.net
￥
257:￥
17/11/20 21:07:09.97 yjl62pX+.net
￥
258:￥
17/11/20 21:07:27.71 yjl62pX+.net
￥
259:￥
17/11/20 21:07:44.86 yjl62pX+.net
￥
260:１３２人目の素数さん
17/11/22 10:59:21.45 1eSvd7ei.net
>>242
検証用のデータがないと誰も追試しないと思う。
>154は模擬データがあったのでやる気になった。
261:￥
17/11/22 12:06:07.41 j3iP9uSb.net
￥
262:￥
17/11/22 12:06:26.73 j3iP9uSb.net
￥
263:￥
17/11/22 12:06:45.80 j3iP9uSb.net
￥
264:￥
17/11/22 12:07:02.42 j3iP9uSb.net
￥
265:￥
17/11/22 12:07:21.82 j3iP9uSb.net
￥
266:￥
17/11/22 12:07:41.50 j3iP9uSb.net
￥
267:￥
17/11/22 12:08:00.58 j3iP9uSb.net
￥
268:￥
17/11/22 12:08:22.66 j3iP9uSb.net
￥
269:￥
17/11/22 12:08:39.11 j3iP9uSb.net
￥
270:￥
17/11/22 12:08:57.12 j3iP9uSb.net
￥
271:１３２人目の素数さん
17/11/24 16:46:01.01 S4G8MO/P.net
関数のオプションの取り得る値を全て表示する方法って知らない？
例えばcor関数のmethodオプションの取り得る値は"pearson", "kendall", "spearman"みたいに。
そういうのは無いのかな？
初めて使う関数のオプションいちいちググるの面倒くさい。
272:￥
17/11/24 18:45:27.78 7RwNGOaZ.net
￥
273:￥
17/11/24 18:45:42.78 7RwNGOaZ.net
￥
274:￥
17/11/24 18:45:58.90 7RwNGOaZ.net
￥
275:￥
17/11/24 18:46:16.24 7RwNGOaZ.net
￥
276:￥
17/11/24 18:46:33.09 7RwNGOaZ.net
￥
277:￥
17/11/24 18:46:49.63 7RwNGOaZ.net
￥
278:￥
17/11/24 18:47:03.93 7RwNGOaZ.net
￥
279:￥
17/11/24 18:47:38.91 7RwNGOaZ.net
￥
280:￥
17/11/24 18:47:57.22 7RwNGOaZ.net
￥
281:￥
17/11/24 18:48:16.46 7RwNGOaZ.net
￥
282:１３２人目の素数さん
17/11/24 19:33:11.50 aKG4T9LS.net
>>265
ヘルプ表示じゃ駄目なのかい？
283:１３２人目の素数さん
17/11/24 20:02:53.10 YIwuPhoP.net
>>265
> args(cor)
function (x, y = NULL, use = "everything", method = c("pearson",
"kendall", "spearman"))
NULL
284:１３２人目の素数さん
17/11/24 22:49:54.91 8tVkv6KQ.net
>>276
前にヘルプで見ようと思ったら書いてなかった覚えが…。マニアックな関数だったんで、そんなもんかな～と諦めてた。やっぱ基本はヘルプか。
>>277
argsは見れる関数と見れない関数があるのであんまり役に立たない印象。
285:￥
17/11/25 00:22:59.73 kb2pRtxG.net
￥
286:￥
17/11/25 00:23:21.36 kb2pRtxG.net
￥
287:￥
17/11/25 00:23:46.66 kb2pRtxG.net
￥
288:￥
17/11/25 00:24:03.98 kb2pRtxG.net
￥
289:￥
17/11/25 00:24:22.25 kb2pRtxG.net
￥
290:￥
17/11/25 00:24:46.55 kb2pRtxG.net
￥
291:￥
17/11/25 00:25:07.98 kb2pRtxG.net
￥
292:￥
17/11/25 00:25:29.17 kb2pRtxG.net
￥
293:￥
17/11/25 00:25:50.80 kb2pRtxG.net
￥
294:￥
17/11/25 00:26:12.77 kb2pRtxG.net
￥
295:１３２人目の素数さん
17/11/25 01:57:31.82 sNF3t9bW.net
>>278
>argsは見れる関数と見れない関数があるのであんまり役に立たない
だったら、例が悪いよ。見られない関数の例をあげないと
> ヘルプで見ようと思ったら書いてなかった
こっちだってわざわざ下請け関数のオプションまでは書かないよ。
ヘルプを見て未解決ならソースを見てくれ。
296:￥
17/11/25 06:17:05.29 kb2pRtxG.net
￥
297:￥
17/11/25 06:17:23.25 kb2pRtxG.net
￥
298:￥
17/11/25 06:17:39.26 kb2pRtxG.net
￥
299:￥
17/11/25 06:17:55.24 kb2pRtxG.net
￥
300:￥
17/11/25 06:18:11.86 kb2pRtxG.net
￥
301:￥
17/11/25 06:18:28.49 kb2pRtxG.net
￥
302:￥
17/11/25 06:18:44.12 kb2pRtxG.net
￥
303:￥
17/11/25 06:19:01.27 kb2pRtxG.net
￥
304:￥
17/11/25 06:19:19.56 kb2pRtxG.net
￥
305:￥
17/11/25 06:19:36.14 kb2pRtxG.net
￥
306:１３２人目の素数さん
17/11/26 17:41:27.37 ahD9OIUm.net
rm(list=ls())
で R のメモリーをクリア
307:１３２人目の素数さん
17/11/28 13:00:14.06 kbsEO2Dt.net
先週長々と質問させていただいたものです。その後、試行錯誤して行った結果十解決し
308:ました。 以前は「summary」としか入れていなかったのですが、「dd$summary」と入力したところ 全体のデータが出力されました。それ以後「dd$summary$dir.area」等で出せるようになりました。 大変初期的な場面で躓いておりましたが、現在何とか進んでおります。 色々とお答えくださった方々、とても参考になりました。本当にありがとうございました。
309:１３２人目の素数さん
17/12/05 15:24:41.67 2L5k1k+Np
semを使用中、specifyModel()でパスを入力後に
NOTE: it is generally simpler to use specifyEquations() or cfa()
see ?specifyEquations
NOTE: adding 3 variances to the model
という文章が出てそのままsem.X <- sem(modelA.X, X, N=n)
と入力すると
solve.default(C[ind, ind]) でエラー:
Lapack routine dgesv: システムは正確に特異です: U[2,2] = 0
追加情報: 警告メッセージ:
1: sem.default(ram, S = S, N = N, raw = raw, data = data, pattern.number = pattern.number, で:
S is not positive-definite: expect problems
2: sem.default(ram, S = S, N = N, raw = raw, data = data, pattern.number = pattern.number, で:
The following variables have no variance or error-variance parameter (double-headed arrow):
x1, x2, x3, x4, x5, x6
The model is almost surely misspecified; check also for missing covariances.
と出るのですがどうすれば警告等が出ずに分析が行えるのでしょうか？
310:１３２人目の素数さん
17/12/10 13:11:49.22 cG4z9akW.net
N(=100)回コインをなげてｎ(=5回)以上続けて表がでる確率。
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)
}
mean(replicate(10^5,seqn()))
> mean(replicate(10^5,seqn()))
[1] 0.81085
案外、高い確率になった。
311:１３２人目の素数さん
17/12/10 22:28:22.98 WuWuRvn7.net
ポアソンでやるやつだっけ
312:１３２人目の素数さん
17/12/11 07:47:07.43 EURogHTh.net
pooledVariance <- function(...) {
args = list(...)
n.args=length(args)
ss2=0
df=0
for(i in 1:n.args){
ss2 = ss2 + var(args[[i]])*(length(args[[i]])-1)
df = df + (length(args[[i]])-1)
}
ss2/df
}
effectsize <- function(y1,y2){
diff=mean(y1)-mean(y2)
var=(var(x1)*(length(x1)-1)+ var(x2)*(length(x2)-1))/(length(c(y1,y2))-2)
sd=sqrt(var)
diff/sd
}
library(effsize)
cohen.d()
313:１３２人目の素数さん
17/12/11 07:50:22.91 EURogHTh.net
>>304
単なる二項分布。
コインが５回続けて表がでたら、0.5^5 <0.05なのでイカサマコインといわれちゃいそうなんだが、
１００回やってみると案外、５回表が続くので確率を計算しようと思ったが、解析的にできる頭がないので
シミュレーションしてみた。
314:１３２人目の素数さん
17/12/11 07:55:42.42 EURogHTh.net
>303は
> rbinom(100,1,0.5)
[1] 0 1 0 0 1 0 0 0 1 1 1 0 1 0 1 1 0 1 1 0 1 1 0 0 0 1 0 1 1 1 1 0
[33] 0 1 0 0 0 0 0 0 1 0 1 1 0 0 0 1 0 0 1 0 1 1 1 1 1 1 1 1 1 1 1 0
[65] 0 0 1 1 0 0 1 0 1 0 0 0 0 0 0 1 0 0 1 1 1 1 1 0 1 1 0 0 1 1 0 0
[97] 0 0 1 1
で５回以上１が連続するときTRUEを返す関数なのだが
もっと簡単にやれないかなぁとは思っている。
rep(1,5)
315:%in% rbinom(100,1,0.5)は １個ずつ評価されてTRUEが５個返ってくるだけ。 文字列にしてgrepを使うとなんとかなりそうな気がしないでもないのだけど....
316:１３２人目の素数さん
17/12/22 17:21:16.51 U1Bq8bdb.net
すまん、教えてほしいだけど
分析するために初めてＲをインストールしようと思って、このスレのあるように公式サイト行ったら、esetが「JS/Redirector.NAV トロイの木馬」を検知したんだが…；
どうしたらいいだ…
317:１３２人目の素数さん
17/12/22 20:09:03.17 ol6XbcBE.net
そっちじゃなくてCRAN行け
318:１３２人目の素数さん
17/12/22 20:42:31.02 U1Bq8bdb.net
ホームページダメとかわけわかんねぇ・・・とんでもねぇな
ありがとう
319:１３２人目の素数さん
17/12/22 20:50:20.52 uzZ0hejh.net
>>308
エロサイトにアクセスしてないｗ？
320:１３２人目の素数さん
17/12/22 21:27:01.31 RNoKRtjI.net
>>310
CRANはここ
URLﾘﾝｸ(cran.r-project.org)
プロジェクトの方はESETの誤検知っぽいんだよな
321:１３２人目の素数さん
17/12/23 09:28:44.66 n0SBd+bp.net
>>311
してねぇよｗｗ
totalvirusで調べると一件引っかかるし、Redirector検知だから、しょうがないね
322:１３２人目の素数さん
17/12/24 20:42:51.55 CT/NKMd7.net
col=rgb(runif(1),runif(1),runif(1),runif(1))で色指定すると
走らせるたびに色がちがっておもしろい。
hist(rnorm(100),col=rgb(runif(1),runif(1),runif(1),runif(1)))
323:１３２人目の素数さん
17/12/31 17:23:19.35 14tdpK/Y.net
stanやJAGSのコードでgamma関数を使おうとして
y = gamma(x)
と、書いたらエラーになった。
stanだと y=tgamma(x)、JAGSだとy=exp(loggam(x))で動作した。
324:１３２人目の素数さん
18/01/02 08:46:50.95 qdmBZ37O.net
ある大学の入学者男女の比率は１であるという帰無仮説を検定する課題が花子と太郎に課された。
花子は５０人を調査できたら終了として入学者を５０人をみつけて１８人が女子であるという結果を得た。
帰無仮説のもとで
５０人中１８人が女子である確率は 0.01603475
これ以下になるのは５０人中０～１８人と３２～５０人が女子の場合なので
両側検定して
> sum(dbinom(c(0:18,32:50),50,0.5))
[1] 0.06490865
> binom.test(18,50,0.5)$p.value
[1] 0.06490865
で帰無仮説は棄却できないと結論した。
URLﾘﾝｸ(i.imgur.com)
一方、本番と十八番が好きな太郎は一人ずつ調べて１８人めの女子がみつかったところで調査を終えることにした。
１８人めがみつかったのは花子と同じく５０人めであった。
帰無仮説のもとで
１８人がみつかるのが５０人めである確率は0.005772512
これ以下になるのは２３人以下５０人以上番めで女子１８人めがみつかった場合なので
両側検定して
pnb=dnbinom(0:999,18,0.5)
> 1 - sum(pnb[-which(pnb<=dnbinom(50-18,18,0.5))]) # < 0.05
[1] 0.02750309
URLﾘﾝｸ(i.imgur.com)
で帰無仮説は棄却される。
どちらの検定が正しいか、どちらも正しくないか？
検定する意図によってｐ値が変わるのは頻度主義統計の欠陥といえるか？
花子の横軸は裏口入学者数、太郎の横軸はサンプル数なので
サンプルでの裏口入学率を横軸にして９５％信頼区間を示す。
花子の検定での信頼区間は0.36～0.72で18/50を含む、p=0.06491
URLﾘﾝｸ(i.imgur.com)
太郎の検定での信頼区間は0.375～0.72で18/50を含まない、p= 0.0275
URLﾘﾝｸ(i.imgur.com)
主観である、検定の中止の基準の差でｐ値や信頼区間が変化するのは変だという批判である。
325:１３２人目の素数さん
18/01/03 10:35:48.54 YJfyxrv+.net
（訂正）
ある大学の入学者男女の比率は１であるという帰無仮説を検定する課題が花子と太郎に課された。
花子は５０人を調査できたら終了として入学者を５０人をみつけて１８人が女子であるという結果を得た。
帰無仮説のもとで
５０人中１８人が女子である確率は 0.01603475
これ以下になるのは５０人中０～１８人と３２～５０人が女子の場合なので
両側検定して
> sum(dbinom(c(0:18,32:50),50,0.5))
[1] 0.06490865
> binom.test(18,50,0.5)$p.value
[1] 0.06490865
で帰無仮説は棄却できないと結論した。
URLﾘﾝｸ(i.imgur.com)
一方、十八という数字が好きな太郎は一人ずつ調べて１８人めの女子がみつかったところで調査を終えることにした。
１８人めがみつかったのは花子と同じく５０人めであった。
帰無仮説のもとで
１８人がみつかるのが５０人めである確率は0.005772512
これ以下になるのは２３人以下５０人以上番めで女子１８人めがみつかった場合なので
両側検定して
pnb=dnbinom(0:999,18,0.5)
> 1 - sum(pnb[-which(pnb<=dnbinom(50-18,18,0.5))]) # < 0.05
[1] 0.02750309
URLﾘﾝｸ(i.imgur.com)
で帰無仮説は棄却される。
どちらの検定が正しいか、どちらも正しくないか？
検定する意図によってｐ値が変わるのは頻度主義統計の欠陥といえるか？
花子の横軸は女子数、太郎の横軸はサンプル数なので
サンプルでの女子の割合を横軸にして９５％信頼区間を示す。
花子の検定での信頼区間は0.36～0.72で18/50を含む、p=0.06491
URLﾘﾝｸ(i.imgur.com)
太郎の検定での信頼区間は0.375～0.72で18/50を含まない、p= 0.0275
URLﾘﾝｸ(i.imgur.com)
主観である、検定の中止の基準の差でｐ値や信頼区間が変化するのは変だという批判である。
326:１３２人目の素数さん
18/01/15 16:30:41.76 wJofbCL/.net
kainokousiki<-function(a,b,c){return (-b+sqrt(b^2-4*a*c))/(2*a)} #解の公式
kainokousiki(1,-5,6)
でrunすると3じゃなくて6を返すんだけど、どこが間違ってる？
327:１３２人目の素数さん
18/01/15 16:35:00.20 wJofbCL/.net
自己解決
かっこが足りなかった
328:１３２人目の素数さん
18/01/19 12:12:23.00 k3h1PrrK.net
Pythonのスレはないのか
329:￥
18/01/21 20:11:34.35 oUqQkvBY.net
￥
330:￥
18/01/21 20:11:56.34 oUqQkvBY.net
￥
331:￥
18/01/21 20:12:22.56 oUqQkvBY.net
￥
332:￥
18/01/21 20:12:47.97 oUqQkvBY.net
￥
333:￥
18/01/21 20:13:05.14 oUqQkvBY.net
￥
334:￥
18/01/21 20:13:34.57 oUqQkvBY.net
￥
335:￥
18/01/21 20:14:02.66 oUqQkvBY.net
￥
336:￥
18/01/21 20:14:16.99 oUqQkvBY.net
￥
337:￥
18/01/21 20:14:35.47 oUqQkvBY.net
￥
338:￥
18/01/21 20:15:07.01 oUqQkvBY.net
￥
339:１３２人目の素数さん
18/02/05 19:08:34.64 AIJT+apj.net
機械学習をきっかけにPythonに逆転された感じだね
340:１３２人目の素数さん
18/02/06 18:35:39.90 tAZA/Fp/.net
『Rを使った～』だとPythonじゃないのかよって思うよね
341:１３２人目の素数さん
18/02/06 21:37:32.10 tUqX17n9.net
ずっとRだけでPython触ったこと無いけど、覚え直す価値ある？
環境構築からもう面倒なイメージ
342:１３２人目の素数さん
18/02/06 22:43:23.82 qQMNyZjW.net
Python自体は:と直後のインデントさえ気を付ければかなり簡単
3系は数が全て小数扱いなので楽
Anacondaというパッケージでインストールすれば、今流行りのJupyter Notebookという開発環境で対話的にコーディングできる（Rも使える）
343:１３２人目の素数さん
18/02/07 20:40:49.48 mgaw9oVw.net
アナコンダてのがRぽくできるのね、ありがとう
dplyrやggplot2みたいに素人でも簡単便利だといいんだけど
Pythonはオブジェクト志向ぽいしすぐ諦めそう
344:１３２人目の素数さん
18/02/12 16:19:11.21 NSJ4iUa4.net
ブラウザ環境なくなっちゃったの？
345:１３２人目の素数さん
18/02/12 17:06:56.33 hMO/kWPg.net
>>336
誤爆？
何をブラウズする環境が無くなったの？
346:１３２人目の素数さん
18/02/13 18:21:54.45 1hic99Cx.net
ブラウザでプログラミングする環境
347:１３２人目の素数さん
18/02/13 21:09:53.07 EFzybxrC.net
>>338
RStudio Serverの話？
348:１３２人目の素数さん
18/02/13 22:27:26.57 Qyl0hNg7.net
RStudio Severあるで。Dockerで使うのがいいんじゃないかな？！
349:１３２人目の素数さん
18/02/25 18:12:28.68 aD34K55o.net
統計学とウェブ解析を交えて実践的な勉強と練習を
したいのですが、おすすめな書籍やサイトはありますか。
実際に解析ツールや分析ツールを用いて
自分で分析解析してから
解答を見て適切な手順や方法、考察を
解説
350:してくれるものが良いです。 統計学は統計検定2級の知識はありますが ウェブ解析はテキスト読んだだけです。
351:１３２人目の素数さん
18/02/26 12:43:55.13 U8DvzUE2.net
ウェブ解析とは具体的になんですか
352:１３２人目の素数さん
18/02/27 00:54:41.17 hS0OJ3qQ.net
R studioは日本語コメント書く度にIMEが無効になったりカーソルがずれたり黒文字の予測変換が黒背景と重なって見えなくなったりと散々だわ
353:１３２人目の素数さん
18/02/27 01:06:58.09 GkJoHhTZ.net
Windows環境で使うから
354:１３２人目の素数さん
18/02/27 01:16:41.39 hS0OJ3qQ.net
会社規定なんでしかたない
UNIX環境使えるのが羨ましい
355:１３２人目の素数さん
18/02/27 01:25:59.40 34tpClZB.net
Docker for WindowsでRStudioサーバー動かせば？
356:１３２人目の素数さん
18/02/27 12:19:33.31 ZcuB0rsn.net
なぜPythonスレッドがない？
357:１３２人目の素数さん
18/02/27 14:05:09.17 O+8uJ5V+.net
板名読める？
358:１３２人目の素数さん
18/03/04 20:14:10.91 nWbPE6nD.net
>>343
IME無効になるのは俺だけじゃなかったと知ってほっとした。
359:１３２人目の素数さん
18/03/04 21:28:34.58 R7ZPBSuG.net
ストアアプリも同じ症状でるからRStudio固有の問題でなくWindows環境の不治の病だと思ってる
360:１３２人目の素数さん
18/03/26 20:57:59.83 WAyucZm1.net
Run any R code you like. There are over three thousand R packages preloaded.
URLﾘﾝｸ(rdrr.io)
361:１３２人目の素数さん
18/03/26 21:11:25.16 WAyucZm1.net
>>351
日本語あると動作しないが、
グラフを描いてくれるのはうれしい。
URLﾘﾝｸ(imagizer.imageshack.com)
362:１３２人目の素数さん
18/04/04 07:52:09.98 PZp+DZN4.net
Fisher test検定時に
p<2.2e-16
と表示されるんですが、これより小さい値の指数桁数を正確に表記する方法教えて下さい。
例えば5.8e-35となるようにです。
363:１３２人目の素数さん
18/04/04 08:05:49.72 lBVECsOD.net
返値の中にあるp値を参照しなされ
364:１３２人目の素数さん
18/04/04 10:02:52.31 CRvlhZKw.net
fisher.test関数の返り値はlist型で、その中にp.valueという名前でp値が格納されているから$演算子を使って直接参照するか、broom::tidy関数に返り値を渡してdata.frame形式で出力してやれば見れる
365:１３２人目の素数さん
18/04/04 12:43:50.08 LbKgW3kd.net
>>353
>>355が言う直接的な参照
> fisher.test(matrix(c(1,120,130,2),2))$p.value
[1] 1.691912e-69
366:１３２人目の素数さん
18/04/04 13:06:17.36 PZp+DZN4.net
352です。よく分かりました！
ありがとうございます！
367:１３２人目の素数さん
18/04/04 14:44:33.05 LbKgW3kd.net
>>357
技術的な助言をしたけど、学術的に言えば、
p < 0.01 は全て p < 0.01 として、具体的なp値を考える意味はないと思うよ。
一部の例外的な研究分野を除いて(e.g., 遺伝統計学)。
368:１３２人目の素数さん
18/04/04 15:18:50.72 PZp+DZN4.net
はい、まさにその例外的な分野で使おうとしてます。ありがとうございます。
369:￥
18/04/07 06:52:38.98 yx+HETs3.net
￥
370:￥
18/04/07 06:53:00.80 yx+HETs3.net
￥
371:￥
18/04/07 06:53:20.77 yx+HETs3.net
￥
372:￥
18/04/07 06:53:40.58 yx+HETs3.net
￥
373:￥
18/04/07 06:54:01.92 yx+HETs3.net
￥
374:￥
18/04/07 06:54:24.32 yx+HETs3.net
￥
375:￥
18/04/07 06:54:47.76 yx+HETs3.net
￥
376:￥
18/04/07 06:55:08.41 yx+HETs3.net
￥
377:￥
18/04/07 06:55:30.46 yx+HETs3.net
￥
378:￥
18/04/07 06:55:54.53 yx+HETs3.net
￥
379:１３２人目の素数さん
18/04/07 10:24:54.55 maRES4NF.net
>>358
多重比較だと意味あるかも
380:１３２人目の素数さん
18/04/19 21:35:35.54 GVMUXyX9.net
Rのガンマ関数はいくつでオーバーフローするかやってみた。
> i=1
> while(gamma(i)!=Inf){
+ i=i+1
+ }
Warning message:
In gamma(i) : value out of range in 'gammafn'
> i
[1] 172
> gamma(172)
[1] Inf
Warning message:
value out of range in 'gammafn'
> gamma(171)
[1] 7.257416e+306
381:１３２人目の素数さん
18/04/29 00:21:51.67 5dW+xNwa.net
matplot()で折れ線グラフ描いたときに、X軸をカテゴリで示したいのですが、
可能でしょうか？
例えばtemp <- c("0時間","8時間","24時間","48時間")として、
matplot()の引数にtempをとるやり方です。
他にもやり方あれば教えてください。
382:１３２人目の素数さん
18/04/29 01:55:59.75 PngyHOQZ.net
>>372
matplot(..., xaxt="n")
axis(1, at=seq(along=temp), lab=temp)
383:１３２人目の素数さん
18/04/30 16:50:51.65 t3vhzyao.net
>>373
遅くなりましたがありがとうございました。
できました！
384:￥
18/04/30 23:55:02.74 y1TqbRSE.net
￥
385:￥
18/04/30 23:55:22.50 y1TqbRSE.net
￥
386:￥
18/04/30 23:55:43.19 y1TqbRSE.net
￥
387:￥
18/04/30 23:56:02.89 y1TqbRSE.net
￥
388:￥
18/04/30 23:56:23.17 y1TqbRSE.net
￥
389:￥
18/04/30 23:56:40.16 y1TqbRSE.net
￥
390:￥
18/04/30 23:56:53.49 y1TqbRSE.net
￥
391:￥
18/04/30 23:57:13.43 y1TqbRSE.net
￥
392:￥
18/04/30 23:57:32.05 y1TqbRSE.net
￥
393:￥
18/04/30 23:57:53.02 y1TqbRSE.net
￥
394:１３２人目の素数さん
18/05/01 18:57:32.34 iUBwAKWd.net
特定の長方形の中に複数の長方形を最小面積で敷き詰める平面充填に関するパッケージってありませんかね
395:１３２人目の素数さん
18/05/01 19:35:44.26 iUBwAKWd.net
やっぱ自分でアルゴリズムを書いてみます
396:１３２人目の素数さん
18/05/02 00:09:39.04 xpXz1bJY.net
Package`deldir'
これどうなの？
397:１３２人目の素数さん
18/05/06 22:15:43.42 BK1CxH7U.net
# jonckheereテストを書いてみた
jonckheere <- function(L,
alternative = c("two.sided", "increasing", "decreasing"),
cat=TRUE){
# L : list of vectors A1,A2,...,Ak, with assumed tendency
How.Many.Greater.Pairs <- function(A,B){ # How many pairs of A[i] > B[j], count as 0.5 when equal,
n.a = length(A)
n.b = length(B)
how.many.greater.pairs = 0
for(i in 1:n.a){
for(j in 1:n.b){
how.many.greater.pairs = how.many.greater.pairs+ifelse(A[i]==B[j],0.5,A[i]>B[j])
}
}
return(how.many.greater.pairs)
}
Sum.of.Greater.Pairs <- function(L){ #L=list(A1,,,,Ak),A1 < A2 < A3,..,< Ak : vector
k = length(L)
comb = combn(1:k,2) # possible combinaition of pairs to compare
n.comb = ncol(comb) # how many combinations
J = 0 # sum of greater pairs
for(i in 1:n.comb){
J = J + How.Many.Greater.Pairs(L[[comb[1,i]]],L[[comb[2,i]]])
}
return(J)
}
J = Sum.of.Greater.Pairs(L)
n = sapply(L,length)
N = sum(n)
EJ = (N^2-sum(n^2))/4
VJ = (N^2*(2*N+3)-sum(n^2*(2*n+3)))/72
Z = (J-EJ)/sqrt(VJ)
alternative = match.arg(alternative)
p.value = switch(alternative, 'two.sided' = 2 * min(pnorm(Z), pnorm(-Z), 0.5),
'increasing' = pnorm(Z),
'decreasing' = pnorm(-Z))
if(cat){
cat( 'p.value = ', p.value,'\n')
cat('alternative hypothesis: ' ,alternative,'\n')
}
invisible(p.value)
}
398:１３２人目の素数さん
18/05/22 20:34:14.53 8iM7y57s.net
代入は「<-」派？「=」派？
399:１３２人目の素数さん
18/05/22 21:16:34.39 iB1pjrmI.net
>>389
実態調査か何か？
<-と=は挙動が違う場合があるので、使い分けていますが、
代入はどっちか�
400:ﾆ問われたら、無論 <- または -> なお、 > 1 -> x これはエラーにならないけど、 > 1 = x 1 = x でエラー: 代入の左辺が不正 (do_set) です これはエラー
401:１３２人目の素数さん
18/05/23 11:17:56.08 OSJ/4EBd.net
>>390
俺は基本=派。
関数の定義は
z.test <- function(x,n=16,sigma=1){
z=sqrt(n)*mean(x)/sigma
2*pnorm(abs(z),lower=FALSE)
}
と書いている。
402:１３２人目の素数さん
18/05/23 11:24:57.67 MGQGuwX9.net
＞1 = x でエラー
当たり前
そんな使い方なんてするかよ
他言語と同じく=一文字の方がすっきりしてイイ
403:１３２人目の素数さん
18/05/24 14:28:38.60 ExPgBsbL.net
こういうのが紛らわしいから、俺は = 推奨。
x <- 1
if(x <- 1) print('YES')
if(x < -1) print('YES')
404:１３２人目の素数さん
18/05/24 14:36:00.11 EQ5K0CF7.net
だよねぇ
<-良くない
405:１３２人目の素数さん
18/05/25 08:05:58.04 ZHt2t+40.net
やったことなかったので関数の初期値設定に<-を使うとどうなるかやってみた。
まず、= の場合
> z.test <- function(x,n=16,sigma=1){
+ z=sqrt(n)*mean(x)/sigma
+ 2*pnorm(abs(z),lower=FALSE)
+ }
> z.test(1:3)
[1] 1.244192e-15
<- で初期値設定すると、エラー
> z0.test <- function(x,n<-16,sigma<-1){
Error: unexpected assignment in "z0.test <- function(x,n<-"
> z=sqrt(n)*mean(x)/sigma
Error in mean(x) : object 'x' not found
> 2*pnorm(abs(z),lower=FALSE)
Error in pnorm(abs(z), lower = FALSE) : object 'z' not found
> }
406:１３２人目の素数さん
18/05/25 08:11:48.02 ZHt2t+40.net
俺は見栄えがいいと思って関数定義には<-を使っているけど
= でも通常に動作する。
> z.test = function(x,n=16,sigma=1){
+ z=sqrt(n)*mean(x)/sigma
+ 2*pnorm(abs(z),lower=FALSE)
+ }
> z.test(1:3)
[1] 1.244192e-15
<- 推奨の人に聞きたいのだけど
<- でないと動作しないってことあるのだろうか？
407:１３２人目の素数さん
18/05/25 09:41:36.92 GCdXSt8a.net
<-は、打つのが
=より面倒
うっとおしい
408:１３２人目の素数さん
18/05/25 20:19:41.68 dXqzteX1.net
RStudioつかってりゃ［Alt+-］で簡単入力
それと引数の指定は代入じゃないと思うんだが感覚が違うよかな？
409:１３２人目の素数さん
18/05/25 21:33:31.44 QNMt6z2O.net
>>397
essだとアンダースコアを入れるを勝手に「<-」になる。
>>395
それは代入ではなく、関数の規定値の設定。
規定値の設定は「=」と決められているので、「<-」はアウト。
ただし、関数を実行するときには「<-」を使うことができる。
> mean(x<-1:10)
[1] 5.5
> x
[1] 1 2 3 4 5 6 7 8 9 10
410:１３２人目の素数さん
18/05/25 22:26:21.26 Cz/+g7h9.net
=に拘り続けて、<-に変える瞬間の君たちがたまらないよ
411:１３２人目の素数さん
18/05/28 16:26:17.30 d/09kgU6.net
>>399
ありがとうございます。
そんな使い方ができたのですね。
こんなのしか知りませんでした。
> f <- function() {
+ x<<-1:10
+ mean (x)
+ }
> f()
[1] 5.5
> x
[1] 1 2 3 4 5 6 7 8 9 10
412:１３２人目の素数さん
18/05/28 20:56:51.59 Osxttqv4.net
一画面しかグラフがないのに
Hit <Return> to see next plot:
と出るの鬱陶しいな、と思っていたら
par(ask=FALSE)
と設定しておけばいいんだな。
413:１３２人目の素数さん
18/05/29 18:08:52.29 ZCWRCH5Y.net
同名のdataがあるのでパッケージを::で指定してもうまくいかなかった。
data(netmeta::parkinson)
でなくて
data(parkinson, package = 'netmeta')
とするんだな。
414:１３２人目の素数さん
18/05/30 06:16:03.29 fLd3NENr.net
hist(rnorm(100),col=rgb(runif(1),runif(1),runif(1),runif(1)))
だと、1色だけど
hist(rnorm(100),col=sample(colours(),(sample(1:10,1))))
なら1～10色で表示される。
415:１３２人目の素数さん
2018/05/31
416:(木) 21:37:03.79 ID:BfVjjX7C.net
417:１３２人目の素数さん
18/05/31 23:16:15.39 grs1zCKo.net
> hist(rnorm(100),col=rgb(runif(1),runif(1),runif(1),runif(1)))
何色に増やしてもいいがな。
何をやりたいのかな？
hist(rnorm(5000),col=apply(matrix(runif(80),4), 2, function(x){rgb(x[1],x[2],x[3],x[4])}))
418:１３２人目の素数さん
18/06/01 21:54:20.67 Ef899k/0.net
:::てどういう時に使うんだろ？
ソース読みたくて
library(BayesFactor)
ttestBF_indepSample
では表示されなかったが、
library(BayesFactor)
BayesFactor:::ttestBF_indepSample
だと出てきた。
419:１３２人目の素数さん
18/06/04 13:04:13.10 3QOAvZa7.net
Stanも覚えなきゃいけないのかよーくそがー
420:１３２人目の素数さん
18/06/04 13:05:22.57 3QOAvZa7.net
::なら名前空間カブールよけよけの術だよね
421:１３２人目の素数さん
18/06/04 14:24:00.43 d7l4LqX9.net
>>408
とりあえずこれ
URLﾘﾝｸ(github.com)
422:１３２人目の素数さん
18/06/05 12:44:14.88 LfwL/Kok.net
>>409
:::と3個なんだよなぁ
423:１３２人目の素数さん
18/06/05 13:06:03.27 CHJNgrVq.net
>>411
グダグダ言っていないで、さっさとヘルプを参照すればいいじゃないか？
?':::'
424:１３２人目の素数さん
18/06/05 19:45:28.11 LfwL/Kok.net
>>412
ありがとう。
stats:::t.test.default
でｔ検定のソースが見られた
425:１３２人目の素数さん
18/06/11 00:14:18.19 K1zGYlRd.net
関数を値で上書きしてしまった場合に元に戻す方法を教えて頂きたく。
例)
rm <- 1
としてしまった場合にrm関数を元に戻したいです。
よろしくお願いします。
426:１３２人目の素数さん
18/06/11 07:21:35.69 +ayC4RGv.net
>>414
直後ならRstudioでctrl+zで元に戻す。
再起動してたらこの方法は無理だろな
427:１３２人目の素数さん
18/06/11 08:24:34.24 2ZHCUqeW.net
>>414
自分で作ったのではない関数の rm だったら rm(rm) で数値ベクトルオブジェクトが消えて、関数の rm が見えるようになる。
428:１３２人目の素数さん
18/06/11 15:01:37.34 mlFyU0v4.net
>>414
> rm <- 1
> ls()
[1] "rm"
> rm(rm)
> ls()
character(0)
> rm
function (..., list = character(), pos = -1, envir = as.environment(pos),
[以下略]
深く考えなくてもrm(rm) で元に戻るけど。
rmが数値ではなくて自作関数だったら、少々ややこしいけど、
base::rm()で大丈夫だろう
429:１３２人目の素数さん
18/06/11 17:04:35.54 Z+okZT62.net
>>414
既に説明されてるけど上書きじゃなくて別オブジェクト扱いになるからrmで消すか明示的にパッケージを::演算子で関数を指定してやればおけ
自作関数だと上書きされちゃうけど
430:１３２人目の素数さん
18/06/11 19:33:14.14 rfEFkwvW.net
rmが自作関数だと戻せないのでは？
> rm = function(x) asin(sqrt(x))
> rm(1)
[1] 1.570796
> rm
function(x) asin(sqrt(x))
> rm = 1
> rm
[1] 1
> rm(rm)
> rm
function (..., list = character(), pos = -1, envir = as.environment(pos),
inherits = FALSE)
431:１３２人目の素数さん
18/06/11 19:51:41.60 Z+okZT62.net
同じ名前空間にあるなら上書きされちゃうよ
.Gloalenv（だったかな？）に自作関数があり、そこに上書きして変数にしてるから消したら消えるだけで戻せない
432:１３２人目の素数さん
18/06/11 19:55:10.50 Z+okZT62.net
ちょっと古いけど、この辺りを読むとなぜそうなるかが分かると思う
（Rの）環境問題について その１。
URLﾘﾝｸ(qiita.com)
433:１３２人目の素数さん
18/06/11 20:04:41.37 mlFyU0v4.net
>>419
戻せます。
> rm <- function(x) asin(sqrt(x))
> rm
function(x) asin(sqrt(x))
> base::rm(rm)
> rm
function (..., list = character(), pos = -1, envir = as.environment(pos),
[以下略]
434:１３２人目の素数さん
18/06/11 20:06:37.76 mlFyU0v4.net
>>420
>.Gloalenv
case sensitiveなRでそれはないのでは？
> help(".GlobalEnv")
435:１３２人目の素数さん
18/06/11 20:19:14.99 Z+okZT62.net
フォローありがとう
436:413
18/06/11 21:00:49.78 K1zGYlRd.net
>>415-423
みなさまご返信ありがとうございました。
例示したものがbase::rmだったので自作関数のことまでは想像だにしませんでした…(がとても良い勉強になりました。)
自環境でも上書きした関数を元に戻すことが出来ることを確認しました。
437:１３２人目の素数さん
18/06/12 07:42:40.88 9RfEtlLW.net
>>422
戻したいのが自作関数
rm = function(x) asin(sqrt(x))
だったら上書きされて無理じゃないの？
ゴミ箱から消去したファイルを復活させるソフトもあるから
PCに詳しければ可能かもしれないけれど
Rの機能としては復活させるのは無理と思う。
438:１３２人目の素数さん
18/06/12 11:00:45.94 cD/mQB+6.net
自作のrm関数はbaseパッケージとは別に定義されるから
base::rm(rm)で自作の方のrmは消えてbase::rmが「復活」します。
439:１３２人目の素数さん
18/06/12 12:07:52.53 cD/mQB+6.net
ああすまん
自作の方を上書きした場合は、そっちの復活は無理ですね。
440:１３２人目の素数さん
18/06/12 14:36:07.63 VF7OcVBc.net
>>426
ちょっと意味がわからない。
自作の関数なら定義を再実行するだけなのに、
消える消えないって何が問題なんだ？
441:１３２人目の素数さん
18/06/13 14:39:27.58 ZXkz+qgi.net
>>425
関数だけでなく既定値でも同じ。
> pi
[1] 3.141593
> pi=2
> pi
[1] 2
> rm(pi)
> pi
[1] 3.14159
442:１３２人目の素数さん
18/06/14 09:27:55.98 mxBGyFKT.net
共同ツール 1
URLﾘﾝｸ(seleck.cc)
URLﾘﾝｸ(trello.com)
ボードのメニュー → Power-Upsから拡張可能 Slack DropBoxなど
Trello Chrome拡張機能 elegant
URLﾘﾝｸ(www.kikakulabo.com)
trelloのオープンソースあり
共同ツール 2
URLﾘﾝｸ(www.google.com)
共同ツール 3
URLﾘﾝｸ(slack.com)
URLﾘﾝｸ(www.dropbox.com)
URLﾘﾝｸ(bitbucket.org)
URLﾘﾝｸ(ja.atlassian.com)
URLﾘﾝｸ(www.sketchapp.com)
URLﾘﾝｸ(photoshopvip.net)
URLﾘﾝｸ(goodpatch.com)
Trello Chrome拡張機能プラグイン集
URLﾘﾝｸ(chrome.google.com)
Slackプラグイン集
URLﾘﾝｸ(slack.com)
Sketchプラグイン集
URLﾘﾝｸ(sketchapp.com)
URLﾘﾝｸ(supernova.studio)
443:１３２人目の素数さん
18/06/16 19:29:24.89 10mr4q2e.net
はっとデシャンｸﾞ
444:１３２人目の素数さん
18/06/17 19:06:33.58 OYjqtCQI.net
>>307
知恵袋に漸化式があったので計算スクリプトを書いてみた。
表がでる確率ｐのコインをN回投げてK回以上表が続く確率。
# N: flips
# K: least sequential head
# p: probability of head
seqNp <- function(N=100,K=5,p=0.5){
if(N==K) return(p^K)
q=1-p
a=numeric(N) # a(n)=P0(n)/p^n , P0(n)=a(n)*p^n
for(i in 1:K) a[i]=q/p^i # P0(i)=q for any i
for(i in K:(N-1)){ # recursive formula
a[i+1]=0 # a[i+1]=q/p*(a[i]+a[i-1]+a[i-2]+...+a[i-(K-1)])
for(j in 0:(K-1)){
a[i+1]=(a[i+1]+a[i-j])
}
a[i+1]=q/p*a[i+1]
}
P0=numeric(N) # P0[n] : probability of ending with 0 head when flipped n times
for(i in 1:N) P0[i]=a[i]*p^i # P0(n)=a(n)*p^n
MP=matrix(rep(NA,N*K),ncol=K)
colnames
445:(MP)=paste0('P',0:(K-1)) MP[,'P0']=P0 MP[1,'P1']=p for(i in (K-2):K) MP[1,i]=0 # MP[k,n] = Pk[n] : probability of ending with k head when flipped n times for(k in 2:K){ # Pk(n+1)=p*P(k-1)(n) for(i in 1:(N-1)) MP[i+1,k]=p*MP[i,k-1] } ret=1-apply(MP,1,sum) ret[N] }
446:１３２人目の素数さん
18/06/18 09:16:00.35 JdABdagZ.net
>>307
この例で、答えとして何が返ってきて欲しいのか？
2で合ってる？
447:１３２人目の素数さん
18/06/18 21:23:58.67 fHaLluDH.net
>>434
1が5回以上出現しているからTRUE（数字なら１）が返り値。
448:１３２人目の素数さん
18/06/19 10:20:04.23 fcA67BpN.net
>>435
run.check <- function (x, n=5) {
# x の中に 1 が n 回以上連続していれば TRUE を返す.
chg <- c(TRUE, diff(x) != 0) # 変化があった場所
chgidx <- c(which(chg), length(x)+1) # 変化があった場所の添え字
run.length <- diff(chgidx) # 0や1の連続している個数
true.length <- run.length[x[chg] == 1] # 1の連続している個数
any(true.length >= n) # 連続している個数が n 以上のrunがあるか?
}
449:１３２人目の素数さん
18/06/19 21:57:17.85 QuKrVrE8.net
>>436
whichとdiffを用いての解法のお返事ありがとうございます。
実行時間を>303のスクリプトと比べてみました。
run.check <- function (x, n=5) {
# x の中に 1 が n 回以上連続していれば TRUE を返す.
chg <- c(TRUE, diff(x) != 0) # 変化があった場所
chgidx <- c(which(chg), length(x)+1) # 変化があった場所の添え字
run.length <- diff(chgidx) # 0や1の連続している個数
true.length <- run.length[x[chg] == 1] # 1の連続している個数
any(true.length >= n) # 連続している個数が n 以上のrunがあるか?
}
# N(=100)回コインをなげてｎ(=5回)以上続けて表がでる確率。
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(10^5,run.check(rbinom(100,1,0.5)))))
user system elapsed
17.56 0.04 17.68
> system.time(mean(replicate(10^5,seqn())))
user system elapsed
9.74 0.07 9.78
後者の方が速いのは1がn個連続したら、TRUEを返して終了するので以後のチェックはしないためだろうと思います。
diff や anyの使い方が勉強になりました。
時間を割いてスクリプトを作成していただいてありがとうございます。<(_ _)>
450:１３２人目の素数さん
18/06/21 10:12:41.96 Ze2kEGUX.net
Rstudioで日本語の入力ができないのですが、解決する方法知ってる方いませんか?
具体的には、IME自体をオンにすることができません。
日本語の表示やコピペは問題なくできます。
環境
Linux - Utubntu
RStudio Version 1.0.143
451:１３２人目の素数さん
18/06/21 10:26:53.84 Fds/4hNR.net
RStudioのQtが古いのでパッチ当てないと日本語入力できないのが現状ですわ
パッチはUbuntu 16.04.LTS用にしかないので他のバージョンならRStudioサーバー使った方が早いかと
なお、パッチはRStudio 日本語とかでグクってみて
452:１３２人目の素数さん
18/06/23 17:38:43.39 aZCZP6wm.net
windows版のRStudioでもたまに日本語入力ができなくなる。
別のソフトでIMEのON/OFFをやってから戻ると直るけど。
453:１３２人目の素数さん
18/06/23 20:18:55.21 4XxgOFIA.net
>>440
うちも同じ。
454:１３２人目の素数さん
18/06/23 21:30:47.03 JcH4CpHc.net
Windowsはストアアプリでも同じ症状でるから持病なんじゃないかと思う
Ubuntuはfcitxのパッチ当てとけばWindowsのような症状は出ないので快適
455:１３２人目の素数さん
18/06/25 13:23:09.19 rmZ6BvQE.net
これをやってみた。
URLﾘﾝｸ(img-cdn.jg.jugem.jp)
スクリプトは
URLﾘﾝｸ(imagizer.imageshack.com)
456:0/t1VcqA.jpg 小数点以下500桁でみると19個素数があった。 > vip=Vectorize(is.prime) > (idx=which (vip(N))) [1] 100 124 150 172 183 202 215 219 255 296 310 314 323 346 400 417 439 441 474 > cat(substring (e,idx[1],idx[1]+9)) 7427466391
457:１３２人目の素数さん
18/06/26 07:08:45.35 Na/Ih9Bj.net
問題
９９人の囚人がいます。彼らの頭に１～１００までのナンバーカードが貼りつけられた帽子をランダムにかぶせます。
他人の帽子は見ることができても、自分の帽子は見ることができません。
帽子の数は全部で１００なので、一つ使われずに余ります。
そのナンバーは囚人達にはわからないようにしておきます。
この状況で、囚人たちに一斉に自分のナンバーを宣言させて、全員が正解だったら釈放するという賭けをします。
囚人たちには帽子をかぶせられる前に相談タイムが設けられています。
どういう戦略を取れば、助かる確率を最も高くできるでしょうか？
URLﾘﾝｸ(000013.blogspot.com)
解答を読んでも数理が理解できないが
確率が0.5になるのはシミュレーションできた。
# URLﾘﾝｸ(000013.blogspot.com)
inversion <- function(x){
n=length(x)
ret=numeric(n)
for(i in 1:(n-1)){
ret[i] = sum(x[i] > x[(i+1):n])
}
sum(ret) # inversion number
}
is.even= function(x) !inversion(x)%%2 # is inverion number even?
prisoner99 <- function(n=100){
indx=sample(1:n,1) # defective number
X=sample((1:n)[-indx])
Y=numeric(n-1)
for (i in 1:(n-1)){ # select as even permutation
x1=X[-i]
x2=(1:n)[!(1:n) %in% x1] # two numbers unseen for i-th prisoner
tmp=X
tmp[i]=x2[1] ; tmp[n]=x2[2]
Y[i]=ifelse(is.even(tmp), x2[1],x2[2])
}
all(X==Y)
}
mean(replicate(1e3,prisoner99()))
458:１３２人目の素数さん
18/07/11 16:20:29.88 zFHa28EV.net
>>438
遅レスだが、Ubuntu 18.04LTSならデフォルトのIBus-mozcで日本語入力ができるそうだ。
詳しくはUbuntu Weekly Recipeの第527回を読んでくだされ
459:１３２人目の素数さん
18/07/22 11:22:38.35 Ott8rTSz.net
覆面算 RED + WHITE = COLOR
どうもＣでやったらうまくいかないのでＲでやってみた。
# RED + WHITE = COLOR
x=c('R','E','D','W','H','I','T','E','C','O','L','O','R')
unique(x)
redwhite <- function(R,E,D,W,H,I,T,C,O,L){
red=10^(2:0)
white=10^(4:0)
color=10^(4:0)
sum(red*c(R,E,D))+sum(white*c(W,H,I,T,E)) - sum(color*c(C,O,L,O,R))
}
x=unique(x) ; x
REDWHITECOLOR <- function(x){
R=x[1]
E=x[2]
D=x[3]
W=x[4]
H=x[5]
I=x[6]
T=x[7]
C=x[8]
O=x[9]
L=x[10]
cat(paste('RED = ',R,E,D,
' WHITE = ',W,H,I,T,E,
' COLOR = ',C,O,L,O,R),'\n')
}
REDWHITECOLOR(unique(x))
library(gtools)
perm=permutations(n=10,r=10,v=0:9)
perm=perm[perm[,1]!=0&perm[,4]!=0&perm[,8]!=0,] # R!=0,W!=0,C!=0
head(perm) ; tail(perm)
n=nrow(perm)
re=numeric(n)
for(i in 1:n){
re[i]=redwhite(perm[i,1],perm[i,2],perm[i,3],perm[i,4],perm[i,5],perm[i,6],perm[i,7],perm[i,8],perm[i,9],perm[i,10])
}
hist(re)
(indx=which(re==0))
REDWHITECOLOR(perm[indx,])
> REDWHITECOLOR(perm[indx,])
RED = 5 8 7 WHITE = 3 9 6 1 8 COLOR = 4
460: 0 2 0 5
461:１３２人目の素数さん
18/07/22 12:50:02.38 Ott8rTSz.net
# AAB
# × CD
# ------
# CCE
# BFAA
# ------
# BEBDE
f <- function(A,B,C,D,E,F){
(A*100+A*10+B)*D==C*100+C*10+E &(A*100+A*10+B)*C==B*1000+F*100+A*10+A &
(C*100+C*10+E)+(B*1000+F*100+A*10+A)*10==B*10000+E*1000+B*100+D*10+E
}
library(gtools)
perm=permutations(n=10,r=6,v=0:9)
n=nrow(perm)
re=numeric(n)
for(i in 1:n){
re[i]=f(perm[i,1],perm[i,2],perm[i,3],perm[i,4],perm[i,5],perm[i,6])
}
indx=which(re==TRUE)
perm[indx,]
462:１３２人目の素数さん
18/07/23 11:09:01.32 SDW/E7TY.net
大変初歩的な質問で恐縮ですが、ご教示いただければ幸いです。
PERT分布を学んでいるのですが、ベータ分布のスケーリングの仕方が分かりません。
最小値a、最頻値ｂ、最大値ｃ として、a, b, c から求めたパラメータをa1、a2 とします。
このとき、
dbeta(x, a1, a2)*(c-a)+a
とすると、当然ながら横軸の値の範囲は0から1で、縦軸の値のみ大きくなってしまいます。
横軸の値の範囲がaからcで、縦軸の値は確率密度のままにするには、どうしたらよいのでしょうか？
463:１３２人目の素数さん
18/07/23 12:15:25.71 llIMc79p.net
え？aだけ右に動かしたいならこれでいいのでは？
dbeta(x-a, a1, a2)
464:１３２人目の素数さん
18/07/23 18:17:20.65 SDW/E7TY.net
>>449
447です。ご教示どうもありがとうございました。
a <- 40 # 最小値
b <- 45 # 最頻値
c <- 70 # 最大値
mu <- (a+4*b+c)/6 # PERT分布の平均値
a1 <- b*(mu-a)/(c-a) # パラメータ a1
a2 <- b*(c-mu)/(c-a) # パラメータ a2
pert1 <- function(x) dbeta(x,a1,a2) # ベータ分布布の計算
まではできたのですが、ご教示いただいた方法で
pert1 <- function(x) dbeta(x-a ,a1,a2)
とすると、ｘ の範囲が0から1で、ｙが0の一様分布のようなグラフになってしまいます。
465:１３２人目の素数さん
18/07/23 20:05:23.59 5MGL8f5L.net
4パラべーたへの変換は下記かな？
URLﾘﾝｸ(en.wikipedia.org)
a1とa2のパラメータ推定って合ってるんだっけ？僕も正直わかんない
やりたいことはこれだと思う
install.packages("mc2d")してからやってね
----------------------------------------------------
rm(list=ls())
library(mc2d)
a <- 40 # 最小値
b <- 45 # 最頻値
c <- 70 # 最大値
pert4<-function(x) dpert(x,min=a,mode=b,max=c,shape=4)
x1<-40:70
plot(x1,pert4(x1))
466:１３２人目の素数さん
18/07/23 21:05:56.52 +HoGXx9d.net
>>451
447です。
450様、まさにこれ、これです。どうもありがとうございました。
おかげさまで大変助かりました。
448様にも、ご助言下さいましたこｔ、厚く御礼申し上げます。
467:１３２人目の素数さん
18/07/27 15:50:12.19 E7dz8jqT.net
lmtestのパッケージをインストールしようとすると不正なマルチバイト文字がありますと出てくるのですがどうすればいいですか？
468:１３２人目の素数さん
18/07/27 17:33:01.36 ufA3ZDTu.net
>>453
Windows 10 （Rの日本語ヘルプはインストールしていない）+ R studioだけど、そのメッセージ出ないな。
469:１３２人目の素数さん
18/07/27 23:29:04.58 klPq25x/.net
ユーザー名が日本語かな
470:１３２人目の素数さん
18/07/28 12:56:05.38 lhV9Awbf.net
職場のSPSS信者たちにRを認めてもらうにはどうしたら良いだろう
「安かろう悪かろう」「俺が知らないものは三流」的な考え方の人が多くて肩身が狭い…
471:１３２人目の素数さん
18/07/28 17:53:14.05 XHIBT/Gx.net
>>456
>>453のようなエラーが出たら「無料のものはやっぱりダメね�
472:v 有償(SPSS) →信頼できる ボランティア→いい加減 と言う思考らしい
473:１３２人目の素数さん
18/07/28 18:33:45.13 l8Bhwsl9.net
>>455
ユーザー名はアルファベット1文字です
studioインストールしてないんですけど、しないとダメですか？
474:１３２人目の素数さん
18/07/28 21:48:13.01 hNHqu0EH.net
他のパッケージはインストールできるの？
475:１３２人目の素数さん
18/07/29 07:46:06.57 u2P99O/7.net
linuxなら環境変数LANGをCにすれば解決する問題だからwindowsでも似たような対処すれば行けそう。
476:１３２人目の素数さん
18/07/29 19:38:42.16 hjP1s1BE.net
これかな？
> Sys.setlocale(locale="C") # locale の変更（USの設定にしないと表示が変になる）
[1] "C"
元に戻すには
> Sys.setlocale(locale="Japanese_Japan.932")
477:１３２人目の素数さん
18/08/02 20:59:18.66 DvtHDake.net
すみません、教えて下さい。
フィッシャーの正確確率検定で、
適合度（一様性）の検定ってできますか？
できるのは独立性の検定のみ？
Rでいろいろ試したみたんですが、
2群ないとエラーになる。。
478:１３２人目の素数さん
18/08/02 21:16:50.41 GWs8Hfcd.net
適合度ならカイ二乗検定でやればいいのでは？
479:１３２人目の素数さん
18/08/02 21:52:11.65 DvtHDake.net
>>463
レス感謝です。
そうなんですが、
期待度数が5未満になってしまった場合は、
カイ二乗検定は向かないんですよね？
480:１３２人目の素数さん
18/08/02 21:54:12.25 82PNyNf0.net
それで他の方法を探しています
481:１３２人目の素数さん
18/08/02 23:02:57.38 clch7kxN.net
>>464
こういうので
r1=5;r2=4;n1=10;n2=12
prop.test(c(r1,r2),c(n1,n2),correct=TRUE)
chisq.test(matrix(c(r1,r2,n1-r1,n2-r2),2,byrow=TRUE),correct=TRUE)
警告がでるという話かな？
Warning message:
In chisq.test(matrix(c(r1, r2, n1 - r1, n2 - r2), 2, byrow = TRUE), :
Chi-squared approximation may be incorrect
482:１３２人目の素数さん
18/08/03 09:04:31.58 VGe+pklE.net
>>464
適合度を見るのに関係ないだろ。
> cnt
[1] 1 0 2 4 3 4 12 6 12 11 21 22 21 37 40 30 59 44 49 47 38 48 36 43 38
[26] 46 33 34 30 21 26 33 17 13 14 12 8 16 7 7 9 6 6 5 3 8 4 4 4 1
[51] 0 2 1 0 0 0 0 2
> prb
[1] 0.0002947255 0.0006464052 0.0012538873 0.0022037831 0.0035720843
[6] 0.0054114364 0.0077413658 0.0105430123 0.0137588816 0.0172972083
[11] 0.0210398698 0.0248524558 0.0285950600 0.0321325300 0.0353432214
[16] 0.0381256526 0.0404028045 0.0421240952 0.0432652773 0.0438266419
[21] 0.0438299769 0.0433147339 0.0423338191 0.0409493611 0.0392287274
[26] 0.0372409836 0.0350539129 0.0327316483 0.0303329196 0.0279098757
[31] 0.0255074215 0.0231629878 0.0209066521 0.0187615256 0.0167443293
[36] 0.0148660914 0.0131329065 0.0115467110 0.0101060386 0.0088067297
[41] 0.0076425765 0.0066058952 0.0056880193 0.0048797146 0.0041715201
[46] 0.0035540195 0.0030180499 0.0025548579 0.0021562071 0.0018144483
[51] 0.0015225575 0.0012741483 0.0010634658 0.0008853653 0.0007352804
[56] 0.0006091856 0.0005035529 0.0004153084
> chisq.test(cnt, p=prb, rescale.p=TRUE, simulat
483:e.p.value=TRUE) Chi-squared test for given probabilities with simulated p-value (based on 2000 replicates) data: cnt X-squared = 65.504, df = NA, p-value = 0.2189 5未満の数があっても関係がない
484:１３２人目の素数さん
18/08/03 10:29:47.68 jRZN4rrt.net
>>467
警告が出るのは独立性検定の場合でしたね。
485:１３２人目の素数さん
18/08/04 02:59:15.87 8dNre62F.net
463です。
カイ二乗検定で、期待度が5未満が多いのは不適切なのかと思い込んでおりましたが、
それは独立性の検定の場合のことであり、
適合度（一様性）の検定の場合は、
気にしなくて良かったのか。。
ならカイ二乗検定でやります。
いろいろ教えて下さり、ありがとうございます。
486:１３２人目の素数さん
18/08/04 19:49:50.94 3IVOmTrE.net
なんでお前らSPSS使わないの？
URLﾘﾝｸ(www.cdleidyba.lt)
487:１３２人目の素数さん
18/08/04 20:37:20.39 kYh4t8ZN.net
>>470
Stanとかに対応してんの？
488:１３２人目の素数さん
18/08/04 22:03:56.61 kTWj7ufe.net
>>470
化石みたいなバージョンやな
489:１３２人目の素数さん
18/08/05 00:56:32.45 biDT5wEJ.net
>>470
いまSPSSを使う理由って何なの？煽りじゃなくマジで
490:１３２人目の素数さん
18/08/05 14:27:54.32 dSJLM/da.net
maple使えよ
491:１３２人目の素数さん
18/08/05 19:10:34.43 mULz5b3r.net
>>473
周囲のspssユーザを見ると「指導教員がspssしか使えないので、
自分もspssしか使えない、
新しい（学術系）ソフトは誰かに丁寧に根気よく指導してみらえないと無理」、
これの再生産。
492:１３２人目の素数さん
18/08/05 19:11:30.48 zaJ/WJYb.net
>>470
研究室に予算がないから
493:１３２人目の素数さん
18/08/05 19:12:42.57 mULz5b3r.net
>>476
予算がなかったら、Rを使え
494:１３２人目の素数さん
18/08/05 19:13:59.36 mULz5b3r.net
>>470と>>476を間違えたorz
495:１３２人目の素数さん
18/08/05 19:33:14.43 mTZcaJ/n.net
予算がないならケーキを食べればいいのに
496:１３２人目の素数さん
18/08/08 13:11:20.17 7/3l/PxD.net
再帰呼び出しに関数名の代わりにRecall使うと処理に要する時間が増えることに気づいた。
フィボナッチ数列の再帰呼び出し
# 関数名での呼び出し
fibo <- function(n){
if(n==1|n==2) return(1)
else fibo(n-1)+fibo(n-2)
}
# Recallを使っての呼び出し
fiboR <- function(n){
if(n==1|n==2) return(1)
else Recall(n-1)+Recall(n-2)
}
> system.time(fibo(30))
user system elapsed
4.27 0.02 4.36
> system.time(fiboR(30))
user system elapsed
6.65 0.00 6.70
497:１３２人目の素数さん
18/08/09 17:47:22.65 l1M8GWWv.net
cumsumを再帰関数で書いてみた。
何度も試行錯誤
# cumsum with recursive call
cumsumR <- function(v){
res=numeric(length(v))
cumsumR_sub <- function(v,res,i){
res[1]=v[1]
if(i==length(v)) return(res)
else{
res[i+1] = res[i] + v[i+1]
Recall(v,res,i+1)
}
}
cumsumR_sub(v,res,1)
}
cumsumR(1:10)
とりえあず、動作した。
> cumsumR(1:10)
[1] 1 3 6 10 15 21 28 36 45 55
498:１３２人目の素数さん
18/08/09 18:31:24.68 l1M8GWWv.net
簡略化できた
cumsumR <- function(v,res=NULL,i=1){
res[1]=v[1]
if(i==length(v)) return(res)
else{
res[i+1] = res[i] + v[i+1]
Recall(v,res,i+1)
}
}
> cumsumR(1:10)
[1] 1 3 6 10 15 21 28 36 45 55
>
499:１３２人目の素数さん
18/08/10 18:41:42.48 Hlm8Oe3x.net
ifelseの動作がどうも納得できないでご教示いただきたいのですが。
これってifelseの仕様でしょうか？バグでしょうか？
> x=2:1
> if(x[1]<=x[2]) x else x[2:1]
[1] 1 2
> ifelse(x[1]<=x[2],x,x[2:1])
[1] 1
500:１３２人目の素数さん
18/08/10 18:55:39.89 vdkgfABT.net
>>483
ifelse() はベクトル演算できる関数。
だから、戻り値の要素数は、test部分の要素数に一致する。
testが1つなので、x[2:1]の1つ目の要素が返される。
501:483
18/08/10 18:58:51.59 vdkgfABT.net
testの要素を2つにすれば、出力も2つになる。
> ifelse(c(x[1] <= x[2], x[1] <= x[2]), x, x[2:1])
[1] 1 2
さらに３つにすると、helpに書いている通り、yesやno部分は再利用される。
> ifelse(c(x[1] <= x[2], x[1] <= x[2], x[1] <= x[2]), x, x[2:1])
[1] 1 2 1
502:１３２人目の素数さん
18/08/10 23:06:17.77 Hlm8Oe3x.net
>>484
ありがとうございました。
バグではなくてそういう仕様だったのですね。
こういう使い方ができるということがわかって勉強になりました。
> x=2:1
> ifelse(c(x[1]!=x[2],x[1]==x[2]),1:2,3:4)
[1] 1 4
503:１３２人目の素数さん
18/08/10 23:17:14.10 VlOUWHrC.net
教わったので早速、
ifelse() はベクトル演算できると、再利用されるの動作確認。
> ifelse(c(TRUE,TRUE,FALSE,FALSE,TRUE,FALSE,FALSE),1:2,11:15)
[1] 1 2 13 14 1 11 12
504:１３２人目の素数さん
18/08/10 23:45:24.37 VlOUWHrC.net
２種類の文字の変換もifelseでできるんだなぁ。
> mat
[,1] [,2] [,3]
[1,] 1 1 1
[2,] 1 1 0
[3,] 1 1 0
[4,] 1 1 0
[5,] 1 1 0
[6,] 1 1 0
> print(apply(mat,2, function(x) ifelse(x==1,'真','偽')),quote=F)
[,1] [,2] [,3]
[1,] 真 真 真
[2,] 真 真 偽
[3,] 真 真 偽
[4,] 真 真 偽
[5,] 真 真 偽
[6,] 真 真 偽
505:１３２人目の素数さん
18/08/11 00:28:09.90 OesSEWnz.net
行列も次元付きのベクトルだからもっと簡略化できた。
> (mat=matrix(sample(0:1,12,rep=T),3,4))
[,1] [,2] [,3] [,4]
[1,] 1 1 0 1
[2,] 0 1 0 0
[3,] 1 0 0 1
> print(ifelse(mat==1,'Right','Wrong'),quote=F)
[,1] [,2] [,3] [,4]
[1,] Right Right Wrong Right
[2,] Wrong Right Wrong Wrong
[3,] Right Wrong Wrong Right
>
506:１３２人目の素数さん
18/08/11 18:04:43.46 05wxJPem.net
RやSPSSって名前の試験まであるのね
507:１３２人目の素数さん
18/08/16 23:26:21.81 3mKQ9SP5+
>>469
適合度検定の場合は、不適の基準は期待度数5未満ではなく、2.5未満では?
つまり一様性の検定なら、全ての期待度数が2.5以上なら良いのでは?
508:１３２人目の素数さん
18/08/22 14:02:15.21 rznk0lAS.net
規約分数にするパッケージが探せなかったので自分でスクリプトを組んでみた。
エラー処置は省略。
reduce_fraction <- function(x,y){
a=x
b=y
r = a%%b # a=bq+r -> a%%b=b%%r
while(r){
a = b
b = r
r = a%%b
}
gcd=b
cat(x/gcd,'/',y/gcd,'\n')
invisible(c(x/gcd,y/gcd))
}
> reduce_fraction(2860,1082900)
11 / 4165
509:１３２人目の素数さん
18/08/23 01:12:39.77 0Wz4IoKN.net
rvest?でログイン必要なとこをスクレイピングする時って
login_page <- html_session("URLﾘﾝｸ(xxxxx"))
login_form <- html_form(login_page)[[1]] %>%
set_values(AAA="xxxx", BBB="xxxx")
というので行けるはずなんですが
最終行のAAAとBBBってそれぞれhtmlタグのname=""のとこからとってくるんですよね?
これにここがname="user[email]"みたいに[]という記号はいってたらどうすればいいでしょうか?
510:１３２人目の素数さん
18/08/25 17:48:59.21 MdsDkupV.net
>>493
よろしくおねがいいしまーーーーーーーす
511:１３２人目の素数さん
18/08/26 22:01:57.06 2bj/HrEX.net
>>494
質問の意味がいまいち分からないから、誰も助言できないのでは？
> x <- 'user[email]'
> x
[1] "user[email]"
rvestパッケージを使ったことがないけど、
記号が入っていたらなぜ問題があるのかよく分からない。
512:１３２人目の素数さん
18/08/26 22:53:43.10 WnqFRbi9.net
というかどこのページのスクレイピングしたいか晒して
仮でもいいから
513:１３２人目の素数さん
18/08/27 06:16:23.20 x8E4FG0O.net
>>495 ちがいます AAAがどれかよくみてみてください ””はありません
515:１３２人目の素数さん
18/08/27 06:19:43.00 x8E4FG0O.net
ちょっと直接はることはできないのでタグだけはります。
<input class="form-control" autofocus="autofocus" placeholder="Email address" required="required" type="text" value="" name="user[email]" id="user_email">>>496
516:１３２人目の素数さん
18/08/27 13:25:39.67 SdOxff6m.net
>>497
もっと分かるように説明しないと、追試できる情報も提供しないし。
もしかして変数名などに[]が入っている場合 にどうしたらよいかってこと？
> `user[email]` <- 1:5
> `user[email]`
[1] 1 2 3 4 5
517:１３２人目の素数さん
18/08/27 15:43:37.46 x8E4FG0O.net
>>499
多分行けました
超感謝!!
518:１３２人目の素数さん
18/08/27 20:20:32.13 5a+trmgU.net
ある問題のシミュレーションしようと思って
問題文の記号のまま
q <- function(x)
....
とやって
q(100)と入力すると
Rが終了することに気づいた。
519:１３２人目の素数さん
18/08/27 22:56:54.68 a+8R0Tu1.net
base::q()
を実行すると爆弾出るね
520:１３２人目の素数さん
18/08/29 22:45:10.53 BcwFyR33.net
ちょっとした疑問です。
空ベクトルの検出って長さ０以外に検出方法ってあるでしょうか？
> x=c(1,2)
> x=x[-1]
> x=x[-1]
> length (x)
[1] 0
> x
numeric(0)
> x==numeric (0)
logical(0)
> is.null(x)
[1] FALSE
> is.na(x)
logical(0)
> x==NULL
logical(0)
> x==NA
logical(0)
> length(x)==0
[1] TRUE
521:１３２人目の素数さん
18/08/31 22:46:15.16 IWQvY6FL.net
４点の座標を入力するとそれらを結ぶ四面体の体積を求めるスクリプトを書いてみた。
高さはパッケージnleqslvを使った近似計算。
# Calculate tetrahedron volume from cordinates
library(nleqslv)
Tetra <- function(O=c(1/2,sqrt(3)/6,sqrt(2/3)),A=c(0,0,0),B=c(1,0,0),C=c(cos(pi/3),sin(pi/3),0)){
fn <- function(x,O,A,B,C){
AO=A-O
BO=B-O
CO=C-O
HO=x[1]*AO+x[2]*BO+(1-x[1]-x[2])*CO # H on triangle ABC
AB=B-A
AC=C-A
c(HO%*%AB,HO%*%AC) # HO vertial to AB and AC
}
fn1 <- function(x) fn(x,O,A,B,C)
x=nleqslv::nleqslv(c(1/3,1/3),fn1)$'x'
AO=A-O
BO=B-O
CO=C-O
HO=x[1]*AO+x[2]*BO+(1-x[1]-x[2])*CO
h=sqrt(sum(HO^2))
a=sqrt(sum((B-C)^2))
b=sqrt(sum((C-A)^2))
c=sqrt(sum((A-B)^2))
s=(a+b+c)/2
base=sqrt(s*(s-a)*(s-b)*(s-c))
V=1/3*base*h
return(V)
}
初期値は辺の長さ１の正四面体
> options(digits = 16)
> Tetra()
[1] 0.1178511301977579
> sqrt(2)/12
[1] 0.1178511301977579
多分、正常に動作していると思う。
522:１３２人目の素数さん
18/09/01 00:25:33.05 52Ub52jp.net
>>504
行列式det使うと簡単
> po <- c(1/2, sqrt(3)/6, sqrt(2/3))
> pa <- c(0,0,0)
> pb <- c(1,0,0)
> pc <- c(cos(pi/3), sin(pi/3), 0)
> det(cbind(pa-po,pb-po,pc-po))/6
[1] -0.1178511
523:１３２人目の素数さん
18/09/01 01:41:47.56 qG52f2Ee.net
ベクトルの三重積を教わったので、パッケージ pracma の外積crossを使った
tetrahedron <- function(O=c(1/2,sqrt(3)/6,sqrt(2/3)),A=c(0,0,0),B=c(1,0,0),C=c(cos(pi/3),sin(pi/3),0)){
AO=A-O
BO=B-O
CO=C-O
as.numeric(abs(pracma::cross(AO,BO) %*% CO)/6)
}
４行で済んだ。
>>505
ありがとうございました。
パッケージに頼らずに計算できたのですね。
524:１３２人目の素数さん
18/09/03 18:55:17.77 S47YTHgP.net
データを解析する前にさらっと特徴を見たい時、皆さんはどんなコマンドを使っていますか？
私が思いつくのは
summary
boxplot
hist
pairs
です。こんなのも良いよってのがあったら教えてくださいm(_ _)m
※ライブラリの使用有無は問いません
525:１３２人目の素数さん
18/09/03 19:38:36.45 fowZfPON.net
>>507
BESTパッケージのplotPost
histに加えて95%CIを表示してくれる。
526:１３２人目の素数さん
18/09/03 20:05:14.66 OGj8hrn2.net
>>507
出力をデータフレーム型にしたいときはskimrパッケージ
527:１３２人目の素数さん
18/09/03 20:07:16.33 OGj8hrn2.net
あと、まだ試してないけsummarytoolsパッケージも面白そう
528:１３２人目の素数さん
18/09/09 22:35:22.23 9XY+z1xx.net
>>503
良い方法が見つからなかったので !length(x)で空白ベクトル判定とした。
文字列を逆に並べる再帰呼び出しスクリプト
> reverse <- function(x){
+ if(!length(x)) return(NULL)
+ c(Recall(x[-1]),x[1])
+ }
> cat(reverse(LETTERS[1:26]))
Z Y X W V U T S R Q P O N M L K J I H G F E D C B A
529:１３２人目の素数さん
18/09/09 23:21:54.06 1udV8jVZ.net
>>511
なぜ if (length(x) == 0) と書かないのか？
可読性って知ってるか？
530:１３２人目の素数さん
18/09/10 16:54:23.64 89eIPezd.net
>>512
文脈読めば
!length(x)で空白ベクトル
の流れ
531:１３２人目の素数さん
18/09/10 18:22:05.74 gLTPxqtw.net
>>513
length は数値を返す関数なのに、なぜわざわざ論理演算をするんだ？
プログラムは、ここの議論を知らない人が読む可能性のほうが高いんだぞ。
532:１３２人目の素数さん
18/09/10 19:12:52.95 ZYY4OYkH.net
>>514
正規分布が一様分布より大きい期待値の算出に
mean(rnorm(100)>runif(100))
って書かない？
俺は
sum(ifelse(rnorm(100)>runif(100),1,0))/100
って書いたりしないけど。
533:１３２人目の素数さん
18/09/10 19:43:25.60 SKQ8XoAj.net
>>514
え、常套手段やん？！
534:１３２人目の素数さん
18/09/10 19:50:35.36 ZYY4OYkH.net
>>516
論理値を数値に置き換えて計算しているから
数値を論理値にしても別に違和感がない。
可読性は慣れの問題。
535:１３２人目の素数さん
18/09/10 19:54:17.32 SKQ8XoAj.net
>>517
アンカー間違えとる？
だから、そんなの常套手段やん？
536:１３２人目の素数さん
18/09/10 21:30:23.37 gLTPxqtw.net
>>515
それとこれとは話が別だ。
logical は TRUE か FALSE の二値しかとらないから、それから数値への変換は自明。
多種の値をとる数値に論理演算をするのはいただけない。
537:１３２人目の素数さん
18/09/10 21:36:20.92 gLTPxqtw.net
>>517
is.empty.vector のような関数があるならそれでよいが、length を使っているのだから、それを勝手にベクトルが空かどうかの判断に使っているのはあなたであって、他人はそのようには考えない、といっているのだ。
自分だけしかこーどを見ないならべつに構わないが、こういうところに晒すのはよくない。
ましてや常套手段などと言って他人に教え込むのはやめてもらいたい。
538:１３２人目の素数さん
18/09/10 21:51:19.40 5QS5/GHY.net
>>520
常套手段というのはmeanの話じゃね？
539:１３２人目の素数さん
18/09/10 21:55:50.32 5QS5/GHY.net
whileの中は１でも２でもよくね？
n=0
while(1){
if(n>10) return(10)
n=n+1
}
n=0
while(2){
if(n>10) return(10)
n=n+1
}
540:１３２人目の素数さん
18/09/10 22:00:57.10 5QS5/GHY.net
>511は
配列を逆順に並べる再帰呼出しのコード。
Rには
541:revという関数がある。 そのコードみてみ！ lengthを真偽判定に使ってますがな。 > base:::rev.default function (x) if (length(x)) x[length(x):1L] else x
542:１３２人目の素数さん
18/09/10 22:18:05.27 DJR6mPp+.net
>>520
>523みた？
if(length(x)!=0)になってないよ。
543:１３２人目の素数さん
18/09/11 00:19:39.45 WPUC3B84.net
>>520
>こういうところ
誰も鵜呑みにしない便所の落書き。。。
544:１３２人目の素数さん
18/09/11 08:48:42.34 QUqp/jpE.net
!を引数が数値のときは0か否かを返す関数と理解すれば
if(!0) print(!1)の結果もサクッとわかる。
545:１３２人目の素数さん
18/09/11 09:08:32.20 QUqp/jpE.net
数値を非零かどうか論理値に変換して処理するかは
関数によるな。
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)
>
546:１３２人目の素数さん
18/09/11 09:12:56.04 QUqp/jpE.net
!ってベクトル対応しているな
> which(!(-1:1))
[1] 2
> !(-2:2)
[1] FALSE FALSE TRUE FALSE FALSE
547:１３２人目の素数さん
18/09/15 08:05:49.52 Vl7XZ52q.net
!をnotではなくis.falseとかis.zero変換すれば( ・∀・)ｲｲ!!
548:１３２人目の素数さん
18/09/15 16:42:00.50 h0gUCZ3r.net
>>529
is.zero だと？
だったら最初から == 0 てよいじゃないか。
549:１３２人目の素数さん
18/09/15 16:45:38.59 h0gUCZ3r.net
>>527
だいたい、if や while は関数じゃないし。
while (1) などと書くのはCに毒されているんじゃねーの？
Rなら while (TRUE) のほうが良いし、repeat というのもある。
550:１３２人目の素数さん
18/09/15 16:45:41.31 Vl7XZ52q.net
>>530
スクリプトを読むときに脳内変換する癖をつけておけば
可読性が悪いと思わずにすむ。
551:１３２人目の素数さん
18/09/15 16:46:13.42 Vl7XZ52q.net
>>520
>523みた？
if(length(x)!=0)になってないよ。
552:１３２人目の素数さん
18/09/15 17:19:35.45 h0gUCZ3r.net
>>533
Rのライブラリ？の書き方なんてあてにならないよ。
そもそも、rev も 509 の reverse も演算の回数を必要最低限にするなら、条件は length(x)) >= 2 または > 1 だ。
まあ if (length(x)) のほうが速かったのかも知れないが、常にそうなるとは限らないのだから、早すぎる最適化は諸悪の根源だ。
553:１３２人目の素数さん
18/09/15 19:30:34.55 Vl7XZ52q.net
可読性の次は速度かよ。
自分の流儀と違う { の使い方でも文句いいそう。
554:１３２人目の素数さん
18/09/15 19:44:27.43 VibLIqgl.net
0,1を1000万個作って
!
!=
>
==
での真偽判定を各々１０回する時間を出してみた。
> 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
555:１３２人目の素数さん
18/09/15 19:57:45.72 VibLIqgl.net
!x と x!=0 で１０回やったが、
> 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
結果は再現された。
556:１３２人目の素数さん
18/09/15 20:06:37.73 VibLIqgl.net
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
御指摘の通り、!の負け
557:１３２人目の素数さん
18/10/04 17:02:48.04 247Ted9r.net
>>445
＞遅レスだが、Ubuntu 18.04LTSならデフォルトのIBus-mozcで日本語入力ができるそうだ。
＞詳しくはUbuntu We
558:ekly Recipeの第527回を読んでくだされ Ubuntu 18.04.1 LTSにアップグレードしたけどできない いままでスクレイピングにつかってrvestのログインコマンドもエラー出るようになったし最悪だ。。。。
559:１３２人目の素数さん
18/10/04 20:04:49.68 W6fWyA5T.net
アップグレードってことは文字入力がFcitxのままになってない？Fcitxならパッチ要だよ。
Voyager 18.04LTS（xubuntu 18.04LTS）ではパッチなしでRStudioに日本語入力できた。
560:１３２人目の素数さん
18/10/26 12:58:59.43 8LsI2NoU.net
昨日、PCを新しくして最新のR入れたんだけど、
作業ディレクトリの固定化できなくね？
前まではショートカットのプロパティで作業フォルダ指定してたんだけど、
新しいR、何度やっても作業ディレクトリがデフォルトのドキュメントに戻りやがる
どうすればいいか、誰がご教示を。
561:１３２人目の素数さん
18/10/26 13:04:50.23 I9GOHEF6.net
>>541
プロパティでコマンドラインオプションを消せばよい。
--cd-ほにゃらら
というのが付いてるでしょ。
562:１３２人目の素数さん
18/10/26 13:13:30.11 8LsI2NoU.net
>>542
あ、なるほど
ありがとうございます、多謝
563:１３２人目の素数さん
18/10/26 23:05:34.11 5yNEfvRx.net
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