# # # ### Creating some demographic variables with 2 REGIONS, and 10 CITIES. # # sample.size.wanted <- 1000 # Insert sample size here, then run the entire script. ################################################################################ # <<< REGION 1 >>> # # Population sizes for each City. n1 <- 165321 n2 <- 281683 n3 <- 261789 n4 <- 200144 n5 <- 185635 N1 <- n1 + n2 + n3 + n4 + n5 sta1 <- n1 sta2 <- sta1 + n2 sta3 <- sta2 + n3 sta4 <- sta3 + n4 sta5 <- sta4 + n5 # Create a Region case/row idenifier variable. region.id <- seq(1:N1) # Create a City identifier varaible and a City name factor. city <- c(rep(1, n1), rep(2, n2), rep(3, n3), rep(4, n4), rep(5, n5)) city.names <- factor(c(rep("El Pasio", n1), rep("Burkley", n2), rep("Ditroit", n3), rep("Bolder", n4), rep("South Beech", n5))) # Create a Gender (or sex) variable. library(car) gender <- factor(rbinom(N1, 1, .53)) gender <- recode(gender, "1 = 'Female'; 0 = 'Male'") # Create an Age variable. age.z <- as.list(0) for (i in 1:5){ n <- c(n1, n2, n3, n4, n5) m <- c(32, 48, 37, 41, 39) sd <- c(13, 29, 18, 22, 20) a <- round(rnorm(n[i], m[i], sd[i])) out <- which(a <= 17) a[out] <- (m[i] + sd[i]) - a[out] upper <- c(105,106,107,108,109,110,110) up <- sample(upper, 1) out <- which(a >= up) a[out] <- a[out] - (m[i] + sd[i]) age.z[[i]] <- a rm(a,out,upper,up) } age <- c(age.z[[1]], age.z[[2]], age.z[[3]], age.z[[4]], age.z[[5]]) rm(age.z, i, n, m, sd) # Create an Education variable. z.age <- (age - mean(age))/sd(age) ed <- .8*z.age - rnorm(N1, 0, .5) z.ed <- (ed - mean(ed))/sd(ed) ed <- (z.ed*2.25) + 11 ed[1:sta1] <- ed[1:sta1] - 3.5 ed[sta1:sta2] <- ed[sta1:sta2] + 1.6 ed[sta2:sta3] <- ed[sta2:sta3] - 2.6 ed[sta3:sta4] <- ed[sta3:sta4] - 1.5 ed[sta4:sta5] <- ed[sta4:sta5] + 1.5 low <- which(age <= 5) ed[low] <- 0 education <- round(ed) rm(ed, low, z.ed) # Create an Income variable. z.edu <- (education - mean(education))/sd(education) z.ic <- rbeta(N1, 3, 18) inc <- .6*z.age + .7*z.edu + .8*city + .95*z.ic income <- (inc * 15000) + 50000 low <- which(income <= 15000) income[low] <- income[low] + 25000 rm(inc, low, z.age, z.edu, z.ic) # Put it all together. region.1 <- data.frame(region.id, city, city.names, gender, age, education, income) rm(region.id, city, city.names, gender, age, education, income, n1, n2, n3, n4, n5, sta1, sta2, sta3, sta4, sta5) detach("package:car") detach("package:MASS") detach("package:nnet") ################################################################################ # <<< REGION 2 >>> # # Population city sizes. n1 <- 331213 n2 <- 119115 n3 <- 262982 n4 <- 208317 n5 <- 95335 N2 <- n1 + n2 + n3 + n4 + n5 sta1 <- n1 sta2 <- sta1 + n2 sta3 <- sta2 + n3 sta4 <- sta3 + n4 sta5 <- sta4 + n5 # Create a Region case/row idenifier variable. region.id <- seq(1:N2) # Create a City identifier varaible and a City factor. city <- c(rep(6, n1), rep(7, n2), rep(8, n3), rep(9, n4), rep(10, n5)) city.names <- factor(c(rep("Whatts", n1), rep("Garvard", n2), rep("Pittsburger", n3), rep("Spingfield", n4), rep("VallStreet", n5))) # Create Gender variables. gender1 <- rbinom(n1, 1, .48) gender2 <- rbinom(n2, 1, .45) gender3 <- rbinom(n3, 1, .49) gender4 <- rbinom(n4, 1, .54) gender5 <- rbinom(n5, 1, .37) # Create an Age variable. age1 <- rnorm(n1) age1 <- round(18 + (age1*5)) low <- which(age1 < 18) age1[low] <- age1[low] + 27 age2 <- rnorm(n2) age2 <- round(30 + (age2*3)) low <- which(age2 < 18) age2[low] <- age2[low] + 30 age3 <- rnorm(n3) age3 <- round(35 + (age3*8)) low <- which(age3 < 18) age3[low] <- age3[low] + 35 age4 <- rnorm(n4) age4 <- round(40 + (age4*9)) low <- which(age4 < 18) age4[low] <- age4[low] + 40 age5 <- rnorm(n5) age5 <- round(32 + (age5*4)) low <- which(age5 < 18) age5[low] <- age5[low] + 28 rm(low) # Create an Education variable. z.age1 <- (age1 - mean(age1))/sd(age1) z.age2 <- (age2 - mean(age2))/sd(age2) z.age3 <- (age3 - mean(age3))/sd(age3) z.age4 <- (age4 - mean(age4))/sd(age4) z.age5 <- (age5 - mean(age5))/sd(age5) zed1 <- (.7*z.age1) + rnorm(n1, 0, 1.0) zed2 <- (.9*z.age2) + rnorm(n2, 0, 1.6) zed3 <- (.7*z.age3) + rnorm(n3, 0, 1.5) zed4 <- (.8*z.age4) + rnorm(n4, 0, 2) zed5 <- (.8*z.age5) + rbinom(n5, 6, .45) ed1 <- round(zed1 + 9) ed2 <- round(zed2 + 18) ed3 <- round(zed3 + 10) ed4 <- round(zed4 + 12) ed5 <- round(zed5 + 14) rm(zed1, zed2, zed3, zed4, zed5) # Create an Income variable. z.edu1 <- (ed1 - mean(ed1))/sd(ed1) z.ic <- rbeta(n1, 3, 18) inc <- .6*z.age1 + .7*z.edu1 + .95*z.ic income1 <- (inc * 2000) + 20000 z.edu2 <- (ed2 - mean(ed2))/sd(ed2) z.ic <- rbeta(n2, 3, 18) inc <- .6*z.age2 + .9*z.edu2 + .95*z.ic income2 <- (inc * 5000) + 50000 z.edu3 <- (ed3 - mean(ed3))/sd(ed3) z.ic <- rbeta(n3, 3, 18) inc <- .8*z.age3 + .6*z.edu3 + .95*z.ic income3 <- (inc * 4000) + 40000 z.edu4 <- (ed4 - mean(ed4))/sd(ed4) z.ic <- rbeta(n4, 8, 40) inc <- .7*z.age4 + .8*z.edu4 + .95*z.ic income4 <- (inc * 6500) + 60000 z.edu5 <- (ed5 - mean(ed5))/sd(ed5) z.ic <- rbeta(n5, 3, 18) inc <- .5*z.age5 + .7*z.edu5 + .95*z.ic income5 <- (inc * 5000) + 95000 rm(inc, z.age1, z.age2, z.age3, z.age4, z.age5, z.edu1, z.edu2, z.edu3, z.edu4, z.edu5, z.ic) # Putting it all together. library(car) gender <- factor(c(gender1, gender2, gender3, gender4, gender5)) gender <- recode(gender, "1 = 'Female'; 0 = 'Male'") rm(gender1, gender2, gender3, gender4, gender5) age <- c(age1, age2, age3, age4, age5) rm(age1, age2, age3, age4, age5) education <- c(ed1, ed2, ed3, ed4, ed5) rm(ed1, ed2, ed3, ed4, ed5) income <- c(income1, income2, income3, income4, income5) rm(income1, income2, income3, income4, income5) region.2 <- data.frame(region.id, city, city.names, gender, age, education, income) rm(region.id, city, city.names, gender, age, education, income, n1, n2, n3, n4, n5, sta1, sta2, sta3, sta4, sta5) detach("package:car") detach("package:MASS") detach("package:nnet") ##### Putting the two Regions together to form the Population. N <- N1 + N2 pop.id <- seq(1:N) region <- factor(c(rep("I", N1), rep("II", N2))) population <- data.frame(pop.id, region, rbind(region.1, region.2)) rm(N, N1, N2, pop.id, region, region.1, region.2); ls() summary(population) ################################################################################ ##### Draw a sample of size n; as specified up at the top of the script. n <- sample.size.wanted s.id <- sample(population[,1], n, replace = F) sample <- population[s.id,-c(3:4)]; rm(n, s.id, sample.size.wanted) names(sample)[1] <- "id" summary(sample) nrow(sample) head(sample) sample <- sample[order(sample\$id),] head(sample) ls() # End of Script.