2017-08-06 81 views
0

感谢您抽出宝贵的时间来阅读。嵌套for循环 - 标记夺回

下面的代码创建的曲线图需要100个样本是5% 和人口(400)的15%之间。

我今天准备这样做,但是,是其它两个部分添加到图表。它 看起来是这样的:

从1-100个样本取100个样本,这些样本在 人口(400)的5%到15%之间。从101-200个人口中抽取100个人口数量在5%到15%之间的样本(800个)。从201-300个人中抽取100个样本,这些样本在5% 和人口的15%(300)之间。

我认为这将需要一个嵌套循环。有没有人有建议作为 如何做到这一点?

谢谢你的时间。基尔斯滕

N <- 400 
pop <- c(1:N) 

lower.bound <- round(x = .05 * N, digits = 0) 
lower.bound ## Smallest possible sample size 

upper.bound <- round(x = .15 * N, digits = 0) 
upper.bound ## Largest possible sample size 

length.ss.interval <- length(c(lower.bound:upper.bound)) 
length.ss.interval ## total possible sample sizes, ranging form lower.bound 
to upper.bound 
sample(x = c(lower.bound:upper.bound), 
     size = 1, 
     prob = c(rep(1/length.ss.interval, length.ss.interval))) 

n.samples <- 100 

dat <- matrix(data = NA, 
       nrow = length(pop), 
       ncol = n.samples + 1) 

dat[,1] <- pop 

for(i in 2:ncol(dat)) { 
    a.sample <- sample(x = pop, 
        size = sample(x = c(lower.bound:upper.bound), 
            size = 1, 
            prob = c(rep(1/length.ss.interval, 
length.ss.interval))), 
        replace = FALSE) 
    dat[,i] <- dat[,1] %in% a.sample 
} 
schnabel.comp <- data.frame(sample = 1:n.samples, 
          n.sampled = apply(X = dat, MARGIN = 2, FUN = 
sum)[2:length(apply(X = dat, MARGIN = 2, FUN = sum))] 
) 
n.prev.sampled <- c(0, rep(NA, n.samples-1)) 
n.prev.sampled 

n.prev.sampled[2] <- sum(ifelse(test = dat[,3] == 1 & dat[,2] == 1, 
           yes = 1, 
           no = 0)) 

for(i in 4:ncol(dat)) { 
    n.prev.sampled[i-1] <- sum(ifelse(test = dat[,i] == 1 & 
rowSums(dat[,2:(i-1)]) > 0, 
            yes = 1, 
            no = 0)) 
} 

schnabel.comp$n.prev.sampled <- n.prev.sampled 
schnabel.comp$n.newly.sampled <- with(schnabel.comp, 
             n.sampled - n.prev.sampled) 
schnabel.comp$cum.sampled <- c(0, 
cumsum(schnabel.comp$n.newly.sampled)[2:n.samples-1]) 
schnabel.comp$numerator <- with(schnabel.comp, 
           n.sampled * cum.sampled) 
schnabel.comp$pop.estimate <- NA 

for(i in 1:length(schnabel.comp$pop.estimate)) { 
    schnabel.comp$pop.estimate[i] <- sum(schnabel.comp$numerator[1:i])/
sum(schnabel.comp$n.prev.sampled[1:i]) 
} 

if (!require("ggplot2")) {install.packages("ggplot2"); require("ggplot2")} 
if (!require("scales")) {install.packages("scales"); require("scales")} 


small.sample.dat <- schnabel.comp 

small.sample <- ggplot(data = small.sample.dat, 
         mapping = aes(x = sample, y = pop.estimate)) + 
    geom_point(size = 2) + 
    geom_line() + 
    geom_hline(yintercept = N, col = "red", lwd = 1) + 
    coord_cartesian(xlim = c(0:100), ylim = c(300:500)) + 
    scale_x_continuous(breaks = pretty_breaks(11)) + 
    scale_y_continuous(breaks = pretty_breaks(11)) + 
    labs(x = "\nSample", y = "Population estimate\n", 
     title = "Sample sizes are between 5% and 15%\nof the population") + 
    theme_bw(base_size = 12) + 
    theme(aspect.ratio = 1) 

我的想法是通过创建嵌套ifelse声明如下:

N.2 <- 800 
N.3 <- 300 
pop.2 <- c(401:N.2) 
pop.3 <- c(801:N) 

lower.bound.2 <- round(x = .05 * N.2, digits = 0) 
upper.bound.2 <- round(x = .15 * N.2, digits = 0) 

lower.bound.3 <- round(x = .05 * N.3, digits = 0) 
upper.bound.3 <- round(x = .15 * N.3, digits = 0) 

也许有些置换...

dat <- imatrix(ifelse(n.samples ,= 100), 
       yes = nrow = length(pop), 
       no = ifelse(n.samples > 100 & > 201), 
       yes = nrow = length(pop.2), 
       no = nrow = length(pop.3), 
       ncol = n.samples + 1) 
+0

是什么意思说“从1-100样本取100个样本是5%的人口(400)的15%之间。” –

+0

如果'N = 400',为什么你要把人口3变成'c(801:N)'? – Onyambu

+0

你可以尝试改说吗?想知道你在做什么是极其困难的 –

回答

0

这是不是你想要做什么?我在下面写的函数mark_recapture有四个参数(样本数,样本的上下限和总体大小),并输出一个矩阵,其中行代表总体中的个体,列代表样本。如果一个人在一个给定的样本中被捕获,它会得到1,否则它会得到0.定义函数后,您可以运行3次,使用3种不同的种群大小来获得3种不同的矩阵。

# define variables 
num_samp <- 100 
lower_sampsize <- 0.05 
upper_sampsize <- 0.15 

# define sampling function that outputs matrix 
mark_recapture <- function (num_samp, pop_size, lower_sampsize, upper_sampsize) { 

    # empty matrix 
    mat <- matrix(0, pop_size, num_samp) 

    # min and max sample size 
    min <- ceiling(lower_sampsize*pop_size) 
    max <- floor(upper_sampsize*pop_size) 

    # vector of random sample sizes between min and max 
    samp_sizes <- sample(min:max, num_samp, replace=TRUE) 

    # draw the samples and fill in the matrix 
    for (i in 1:num_samp) {mat[sample(1:pop_size, samp_sizes[i]),i] <- 1} 

    # return matrix 
    return(mat) 
} 

# do the sampling from the 3 populations 
mat1 <- mark_recapture(num_samp=num_samp, pop_size=400, lower_sampsize=lower_sampsize, upper_sampsize=upper_sampsize) 
mat2 <- mark_recapture(num_samp=num_samp, pop_size=800, lower_sampsize=lower_sampsize, upper_sampsize=upper_sampsize) 
mat3 <- mark_recapture(num_samp=num_samp, pop_size=300, lower_sampsize=lower_sampsize, upper_sampsize=upper_sampsize) 

虽然超出了这个问题的范围,我只提到,有专门的R程序包,以分析和模拟标志重捕的数据,例如,multimark。只需Google“CRAN标记重新获得”即可找到许多选项。我会建议仔细观​​察一下,仔细想想你在这里做什么,因为你可能试图重新发明轮子。