[R] Sampling from a population
Mario Sendra Pina
Mario.Sendra at uv.es
Tue Nov 20 17:44:05 CET 2001
Suppose you have a population of N <- 5 observations, x <- c(43, 28, 7, 61, 39). From that you can draw a maximum of 10 samples without replacement of size n <- 3. (Command choose(N,n) yields 10). For instance the samples I seek are
43, 61, 7
39, 7, 28 ...etc
How can I get R to do that for me, to get an exhaustive list of samples of size n drawn without replacement from a population of size N? The command, sample(x, 3, replace=FALSE), works well for one draw. Is there a package that will handle multiple draws?
# First solution
combix<-function(x,r) {
n<-length(x)
total<-choose(n,r)
a<-matrix(0,total,r)
ind<-c(1:r)
for (i in 1:total){
for (j in 1:r){a[i,j]<-x[ind[j]]}
for (j in r:1){
miro<-j
ind[j]<-ind[j]+1
if (ind[j]<n-r+j+1) break()
}
while(miro<r) {
ind[miro+1]<-ind[miro]+1
miro<-miro+1
}
}
return(a)
}
print(combix(x=1:5,r=3))
x<-c(8,3,1,11,4,7)
print(mean(x))
medias<- apply(combix(x,2),1,mean)
print(mean(medias))
#With-Replacement
combixr<-function(x,r) {
n<-length(x)
total<-n**r
a<-matrix(0,total,r)
ind<-rep(1,r)
for (i in 1:total){
for (j in 1:r){a[i,j]<-x[ind[j]]}
for (j in r:1){
miro<-j
ind[j]<-ind[j]+1
if (ind[j]<n+1) break()
}
while(miro<r) {
ind[miro+1]<-1
miro<-miro+1
}
}
return(a)
}
print(combixr(x=1:5,r=3))
x<-c(8,3,1,11,4,7)
print(mean(x))
medias<- apply(combixr(x,2),1,mean)
mean(medias)
#Other form
darcombi<- function(N,n,caso){
vecpeq<-rep(0,n)
ante<-0
ante1<-0
for(i in 1:n){
N1<- N-(ante1+1)
n1<- n-i
v1<-cumsum(choose(N1:n1,n1))
ante<-min(which(caso<=v1))
vecpeq[i]<-ante+ante1
resta<-0
if (ante>1) resta<-v1[ante-1]
caso<-caso-resta
ante1<-vecpeq[i]
}
return(vecpeq)
}
N<-5
n<-3
for (i in 1:choose(N,n)){
print(darcombi(N=N,n=n,i))
}
x <- c(43, 28, 7, 61, 39)
N<-length(x)
n<-3
for (i in 1:choose(N,n)){
print(x[darcombi(N=N,n=n,i)])
}
-------------- next part --------------
An HTML attachment was scrubbed...
URL: https://stat.ethz.ch/pipermail/r-help/attachments/20011120/5dd9b6ff/attachment.html
More information about the R-help
mailing list