R Code for PPS Sampling

Installing Packages

#install.packages("pps")
#install.packages("survey")
#install.packages("sampling")
#install.packages("writexl")
#install.packages("dplyr")
#install.packages("srvy")

Loading Packages

library(readxl)

## Warning: package 'readxl' was built under R version 4.5.3

library(survey)

## Warning: package 'survey' was built under R version 4.5.3

## Loading required package: grid

## Loading required package: Matrix

## Loading required package: survival

## 
## Attaching package: 'survey'

## The following object is masked from 'package:graphics':
## 
##     dotchart

library(dplyr)

## Warning: package 'dplyr' was built under R version 4.5.3

## 
## Attaching package: 'dplyr'

## The following objects are masked from 'package:stats':
## 
##     filter, lag

## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union

library(sampling)

## Warning: package 'sampling' was built under R version 4.5.3

## 
## Attaching package: 'sampling'

## The following objects are masked from 'package:survival':
## 
##     cluster, strata

library(srvyr)

## Warning: package 'srvyr' was built under R version 4.5.3

## 
## Attaching package: 'srvyr'

## The following object is masked from 'package:stats':
## 
##     filter

library(pps)

Reading Data

base<-readRDS(file = "Marco.rds")

Sampling Frame of Municipalities

municipios=base%>%group_by(COLE_COD_MCPIO_UBICACION,COLE_MCPIO_UBICACION)%>%
  summarise(x=n(),y=mean(PUNT_GLOBAL))

## `summarise()` has regrouped the output.
## ℹ Summaries were computed grouped by COLE_COD_MCPIO_UBICACION and
##   COLE_MCPIO_UBICACION.
## ℹ Output is grouped by COLE_COD_MCPIO_UBICACION.
## ℹ Use `summarise(.groups = "drop_last")` to silence this message.
## ℹ Use `summarise(.by = c(COLE_COD_MCPIO_UBICACION, COLE_MCPIO_UBICACION))` for
##   per-operation grouping (`?dplyr::dplyr_by`) instead.

Selection Probabilities of Municipalities

municipios$pk=municipios$x/sum(municipios$x)

Selecting an Ordered Sample

m=50
set.seed(123)
s<-ppswr(municipios$x,m)
muestra_mun=municipios[s,]
muestra_mun=muestra_mun[order(-muestra_mun$x),]

Creating the “dsgn” Objec

dsgn_mun=svydesign(id=~1,data=muestra_mun,probs=~pk)

Estimating the Global Mean Score

(est=svymean(~y,dsgn_mun,deff=T))

##       mean       SE  DEff
## y 245.7802   2.4033 2.845

salida=function(est,alpha){
  est=as.data.frame(est)
  names(est)[2]="se"
  est$cv=100*(est$se/est[,1])
  est$ic_low=est[,1]-qnorm(1-alpha/2)*est$se 
  est$ic_upp=est[,1]+qnorm(1-alpha/2)*est$se
  return(round(est,2))
}

alpha=0.05

(tabla=salida(est,alpha))

##     mean  se deff   cv ic_low ic_upp
## y 245.78 2.4 2.84 0.98 241.07 250.49