library(TeachingSampling)
library(tidyverse)
library(e1071)
library(fitdistrplus)Técnicas de Amostragem - Atividade Computacional
Departamanto de Estatística e Matemática Aplicada - UFC
1 Base de Dados Lucy
Este conjunto de dados corresponde a uma amostra aleatória do BigLucy. Ele contém algumas variáveis financeiras de 2396 empresas industriais de uma cidade em um determinado ano fiscal.
As variáveis são:
ID: O identificador da empresa. Corresponde a uma sequência alfanumérica (duas letras e três dígitos).
Ubication (Localização): O endereço do escritório principal da empresa na cidade.
Level (Nível): As empresas industriais são discriminadas de acordo com os impostos declarados. Existem empresas pequenas, médias e grandes.
Zone (Zona): A cidade é dividida em zonas geográficas. Uma empresa é classificada em uma zona específica de acordo com seu endereço.
Income (Renda): O montante total dos ganhos (ou lucros) de uma empresa no ano fiscal anterior. É calculado levando em consideração as receitas e ajustando os custos operacionais.
Employees (Funcionários): O número total de pessoas trabalhando para a empresa no ano fiscal anterior.
Taxes (Impostos): O montante total do Imposto de Renda de uma empresa.
SPAM: Indica se a empresa utiliza as opções de Internet e WEBmail para fazer autopropaganda.
A tabela 1 apresenta os 15 primeiros registros da base de dados:
data("Lucy")
attach(Lucy)
df = as.data.frame(Lucy)
head(df, 15) |> knitr::kable(caption = "15 primeiros registros do
dataset Lucy")| ID | Ubication | Level | Zone | Income | Employees | Taxes | SPAM |
|---|---|---|---|---|---|---|---|
| AB001 | c1k1 | Small | A | 281 | 41 | 3.0 | no |
| AB002 | c1k2 | Small | A | 329 | 19 | 4.0 | yes |
| AB003 | c1k3 | Small | A | 405 | 68 | 7.0 | no |
| AB004 | c1k4 | Small | A | 360 | 89 | 5.0 | no |
| AB005 | c1k5 | Small | A | 391 | 91 | 7.0 | yes |
| AB006 | c1k6 | Small | A | 296 | 89 | 3.0 | no |
| AB007 | c1k7 | Small | A | 490 | 22 | 10.5 | yes |
| AB008 | c1k8 | Small | A | 473 | 57 | 10.0 | yes |
| AB009 | c1k9 | Small | A | 350 | 84 | 5.0 | yes |
| AB010 | c1k10 | Small | A | 361 | 25 | 5.0 | no |
| AB011 | c1k11 | Small | A | 374 | 34 | 6.0 | yes |
| AB012 | c1k12 | Small | A | 419 | 20 | 7.0 | no |
| AB013 | c1k13 | Small | A | 402 | 19 | 7.0 | yes |
| AB014 | c1k14 | Small | A | 330 | 23 | 4.0 | yes |
| AB015 | c1k15 | Small | A | 411 | 36 | 7.0 | no |
A análise será focada nas variáveis Level (Nível), Zone (Zona), Income (Renda) e SPAM. Inicialmente, será realizada uma análise exploratória das variáveis selecionadas, em seguida será realizada uma amostragem Bernoulli e tentaremos estimar a renda total da população. Por último, irei sujerir uma estratégia amostral a partir da análise realizada inicialmente.
1.1 1. Análise Descritiva
Dentre as variáveis selecionadas apenas a variável Income (Renda) é uma variável quantitativa (contínua), as variável Level (Nível) é qualitativa ordinal uma vez que ela possui uma ordenação hierárquica (Small, Medium e Big), a variável Zone (Zona) é qualitativa e apresenta cinco categorias (A, B, C, D e E), já a variável e SPAM é qualitativa de característica dicotômica, uma vez que apresenta dois possíveis valores, sendo eles yes ou no.
Para possíveis planos amostrais é importante considerar a quantidade de empresas localizadas em cada uma das cinco zonas (A, B, C, D e E). Das 2396 empresas, 307 estão localizadas na Zona A, 727 estão na Zona B, 974 estão na Zona C, 223 estão na Zona D e 165 estão localizadas na Zona E.
df |>
ggplot(aes(x = Zone)) +
geom_bar(fill = "lightblue3") +
geom_text(stat = 'count', aes(label = ..count..), vjust = -0.5) +
labs(x = "Zona", y = "") +
theme_minimal() +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank()) +
coord_cartesian(ylim = c(0, 1000))
Ademais, quando analisamos a quantidade de empresas divididas pelos seus níveis de tamanho, observamos que das 2396 empresas, 1576 são classificadas como pequenas, 737 são classificadas como médias e 83 empresas são consideradas grandes.
df |>
count(Zone) |>
ggplot(aes(2, n, fill = Zone)) +
geom_col() +
geom_text(aes(label = paste0(round((n/sum(n))*100, 1), "%")),
position = position_stack(vjust = 0.5),
color = "black", size = 3) +
scale_fill_manual(values = c("coral1", "darkseagreen3", "lightblue3",
"lightgoldenrod1", "darkgray")) +
coord_polar("y") +
xlim(c(0.2, 3)) +
theme_minimal() +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank()) +
labs(title = "Porcentagem de Empresas por Zonas")
df |>
ggplot(aes(x = Level)) +
geom_bar(fill = "lightblue3") +
geom_text(stat = 'count', aes(label = ..count..), vjust = -0.5) +
labs(x = "Nível", y = "") +
theme_minimal() +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank()) +
coord_cartesian(ylim = c(0, 1600))
Das 307 empresas presentes na zona A, 97 delas são classificadas como pequenas, 180 como médias e 30 são classificadas como grandes.
Das 727 empresas presentes na zona B, 593 delas são classificadas como pequenas, 121 como médias e 13 são classificadas como grandes.
Das 974 empresas presentes na zona C, 862 delas são classificadas como pequenas, 111 como médias e apenas 1 é classificada como grande.
Das 223 empresas presentes na zona D, 20 delas são classificadas como pequenas, 187 como médias e 6 são classificadas como grandes.
Das 165 empresas presentes na zona E, 4 delas são classificadas como pequenas, 138 como médias e 23 são classificadas como grandes.
df |> count(Level) |> ggplot(aes(2, n, fill = Level)) +
geom_col() +
labs(x = '',
y = '',
fill = 'Nível') +
scale_fill_manual(values = c("coral1", "darkseagreen3", "lightblue3")) +
coord_polar("y") +
xlim(c(0.2, 3)) +
geom_col() +
geom_text(aes(label = paste0(round((n/sum(n))*100, 1), "%")),
position = position_stack(vjust = 0.5),
color = "black", size = 3) +
theme_minimal() +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank()) +
labs(title = "Porcentagem de Empresas por Nível")
df |>
ggplot(aes(y = fct_rev(fct_infreq(Zone)), fill = Level)) +
geom_bar(width = .8) +
labs(x = '',
y = '',
fill = 'Nível') +
scale_fill_manual(values = c("coral1", "darkseagreen3", "lightblue3")) +
geom_text(aes(label = after_stat(count)),
stat = 'count',
position = position_stack(vjust = 0.25),
size = 3) +
theme_minimal() +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank()) +
coord_cartesian(xlim = c(0, 1000))
1.2 SPAM
A variável SPAM nos dá uma ideia de quantas empresas do dataset Lucy investiram em algum tipo de propaganda na internet ou webmail para se autopromoverem. Observamos que, das 2396 empresas, 1459 investiram e 937 não investiram.
df |>
ggplot(aes(x = SPAM)) +
geom_bar(fill = "lightblue3") +
geom_text(stat = 'count', aes(label = ..count..), vjust = -0.5) +
labs(x = "SPAM", y = "") +
theme_minimal() +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank()) +
coord_cartesian(ylim = c(0, 1500))
A seguir, estão as principais medidas descritivas da variável Renda. Obervamos uma renda média nessa amostra de 432,0605 com um desvio padrão de 266,9792 e um coeficiente de varição de 61,79%.
df |>
summarise(
Mean = mean(Income),
"1st Qu." = quantile(Income, probs = 0.25),
Median = median(Income),
"3st Qu." = quantile(Income, probs = 0.75),
Desv.pd = sd(Income),
Min = min(Income),
Max = max(Income),
Kurt = kurtosis(Income),
Skew = skewness(Income),
CV = sprintf("%.2f%%", (sd(Income)/mean(Income))*100)
) |> knitr::kable(caption = "Medidas descritivas da variável
Income (Renda)", digits = 4)| Mean | 1st Qu. | Median | 3st Qu. | Desv.pd | Min | Max | Kurt | Skew | CV |
|---|---|---|---|---|---|---|---|---|---|
| 432.0605 | 230 | 390 | 576 | 266.9792 | 1 | 2510 | 3.7451 | 1.3834 | 61.79% |
A variável Income pode ser modelada por meio de uma distribuição Gamma de parâmetros 2,570683886 (escala) e 0,005949936 (forma):
library(fitdistrplus)
fitgamma = df$Income |> fitdist("gamma", method = "mle");fitgammaFitting of the distribution ' gamma ' by maximum likelihood
Parameters:
estimate Std. Error
shape 2.570683886 0.0644231357
rate 0.005949936 0.0001598176df |>
ggplot(aes(x = Income)) +
geom_histogram(width = .25, fill = "lightblue3",
na.rm = TRUE, aes(y = after_stat(density))) +
stat_function(fun=dgamma,
args=list(shape = 2.570683886,
rate = 0.005949936),
colour="red") +
labs(x = "Renda", y = "") +
theme_minimal() +
theme(
panel.grid.major = element_blank(),
panel.grid.minor = element_blank())
Quando analisamos a distribuição dos dados da variável renda, notamos uma assimetria à direita, o que é comprovado pela medida de assimetria calculada 1,3834.
Entretanto, isso possui limitações para nossa análise, uma vez que não é possível ter uma ideia segmentada por zonas ou níveis de empresas. Consequantemente, analisar apenas a variável renda isoladamente pode trazer resultados desbalanceados e ineficientes, uma vez que ela por sí só não capta a variação de categorias existentes na amostra do dataset Lucy.
Dessa maneira, será realizada uma análise exploratória levando em consideração as outras variáveis selecionadas, para assim, termos uma ideia de como se comporta a renda quando levamos em consideração os tamanhos das empresas e as zonas que elas estão localizadas.
1.3 Renda por Nível da empresa
df |> filter(Level == "Small") |>
summarise(Mean = mean(Income),
"1st Qu." = quantile(Income, probs = 0.25),
Median = median(Income),
"3st Qu." = quantile(Income, probs = 0.75),
Desv.pd = sd(Income),
Min = min(Income),
Max = max(Income),
Kurt = kurtosis(Income),
Skew = skewness(Income),
CV = sprintf("%.2f%%", (sd(Income)/mean(Income))*100)) |>
knitr::kable(caption = 'Medidas descritivas da variável Income (Renda)
das empresas pequenas',
digits = 4)| Mean | 1st Qu. | Median | 3st Qu. | Desv.pd | Min | Max | Kurt | Skew | CV |
|---|---|---|---|---|---|---|---|---|---|
| 281.8274 | 181 | 280 | 385 | 123.053 | 1 | 499 | -1.0817 | -0.0104 | 43.66% |
df |> filter(Level == "Medium") |>
summarise(Mean = mean(Income),
"1st Qu." = quantile(Income, probs = 0.25),
Median = median(Income),
"3st Qu." = quantile(Income, probs = 0.75),
Desv.pd = sd(Income),
Min = min(Income),
Max = max(Income),
Kurt = kurtosis(Income),
Skew = skewness(Income),
CV = sprintf("%.2f%%", (sd(Income)/mean(Income))*100)) |>
knitr::kable(caption = 'Medidas descritivas da variável Income (Renda)
das empresas médias',
digits = 4)| Mean | 1st Qu. | Median | 3st Qu. | Desv.pd | Min | Max | Kurt | Skew | CV |
|---|---|---|---|---|---|---|---|---|---|
| 661.2632 | 560 | 630 | 740 | 126.8144 | 500 | 999 | -0.2083 | 0.8421 | 19.18% |
df |> filter(Level == "Big") |>
summarise(Mean = mean(Income),
"1st Qu." = quantile(Income, probs = 0.25),
Median = median(Income),
"3st Qu." = quantile(Income, probs = 0.75),
Desv.pd = sd(Income),
Min = min(Income),
Max = max(Income),
Kurt = kurtosis(Income),
Skew = skewness(Income),
CV = sprintf("%.2f%%", (sd(Income)/mean(Income))*100)) |>
knitr::kable(caption = 'Medidas descritivas da variável Income (Renda)
das empresas grandes',
digits = 4)| Mean | 1st Qu. | Median | 3st Qu. | Desv.pd | Min | Max | Kurt | Skew | CV |
|---|---|---|---|---|---|---|---|---|---|
| 1249.47 | 1078.5 | 1153 | 1393 | 253.7679 | 1003 | 2510 | 5.8536 | 1.9469 | 20.31% |
df |> ggplot(aes(Level, Income)) +
stat_boxplot(geom = 'errorbar', na.rm = T, width = .25) +
geom_boxplot(na.rm = T) +
theme_test()
Quando comparamos a renda das empresas entre as categorias Small, Medium e Big, observamos que renda médias das empresas pequenas é de 281,8274, das empresas médias é de 661,2632 e das empresas grandes é de 1249,47. Ou seja, a renda das empresas consideradas grandes é mais que o dobro das empresas médias e quatro vezes maior do que das empresas pequenas. Além disso, um coeficiente de variação (CV) de 43,66% das empresas pequenas, contra 19,18% das empresas médias e 20,31% das empresas grandes indica uma maior variação na renda das empresas pequenas. Por meio do boxplot acima é possível perceber visualmente a disperção da renda por nível da empresa. A categoria Big é a única que possui outliers.
Iremos, agora, fazer uma análise de como se comporta a renda quando comparamos em relação à variável Zone
1.4 Renda por Zona
df |> filter(Zone == "A") |>
summarise(Mean = mean(Income),
"1st Qu." = quantile(Income, probs = 0.25),
Median = median(Income),
"3st Qu." = quantile(Income, probs = 0.75),
Desv.pd = sd(Income),
Min = min(Income),
Max = max(Income),
Kurt = kurtosis(Income),
Skew = skewness(Income),
CV = sprintf("%.2f%%", (sd(Income)/mean(Income))*100)) |>
knitr::kable(caption = 'Medidas descritivas da variável Income (Renda)
das empresas da zona A')| Mean | 1st Qu. | Median | 3st Qu. | Desv.pd | Min | Max | Kurt | Skew | CV |
|---|---|---|---|---|---|---|---|---|---|
| 652.2834 | 450 | 610 | 768.5 | 281.4491 | 210 | 2510 | 6.244997 | 1.717336 | 43.15% |
df |> filter(Zone == "B") |>
summarise(Mean = mean(Income),
"1st Qu." = quantile(Income, probs = 0.25),
Median = median(Income),
"3st Qu." = quantile(Income, probs = 0.75),
Desv.pd = sd(Income),
Min = min(Income),
Max = max(Income),
Kurt = kurtosis(Income),
Skew = skewness(Income),
CV = sprintf("%.2f%%", (sd(Income)/mean(Income))*100)) |>
knitr::kable(caption = 'Medidas descritivas da variável Income (Renda)
das empresas da zona B')| Mean | 1st Qu. | Median | 3st Qu. | Desv.pd | Min | Max | Kurt | Skew | CV |
|---|---|---|---|---|---|---|---|---|---|
| 320.7469 | 188 | 262 | 338.5 | 227.2832 | 1 | 1616 | 5.150401 | 1.986357 | 70.86% |
df |> filter(Zone == "C") |>
summarise(Mean = mean(Income),
"1st Qu." = quantile(Income, probs = 0.25),
Median = median(Income),
"3st Qu." = quantile(Income, probs = 0.75),
Desv.pd = sd(Income),
Min = min(Income),
Max = max(Income),
Kurt = kurtosis(Income),
Skew = skewness(Income),
CV = sprintf("%.2f%%", (sd(Income)/mean(Income))*100)) |>
knitr::kable(caption = 'Medidas descritivas da variável Income (Renda)
das empresas da zona C')| Mean | 1st Qu. | Median | 3st Qu. | Desv.pd | Min | Max | Kurt | Skew | CV |
|---|---|---|---|---|---|---|---|---|---|
| 331.0195 | 189.25 | 362.5 | 450 | 153.5773 | 20 | 1360 | 1.030268 | 0.2515107 | 46.40% |
df |> filter(Zone == "D") |>
summarise(Mean = mean(Income),
"1st Qu." = quantile(Income, probs = 0.25),
Median = median(Income),
"3st Qu." = quantile(Income, probs = 0.75),
Desv.pd = sd(Income),
Min = min(Income),
Max = max(Income),
Kurt = kurtosis(Income),
Skew = skewness(Income),
CV = sprintf("%.2f%%", (sd(Income)/mean(Income))*100)) |>
knitr::kable(caption = 'Medidas descritivas da variável Income (Renda)
das empresas da zona D')| Mean | 1st Qu. | Median | 3st Qu. | Desv.pd | Min | Max | Kurt | Skew | CV |
|---|---|---|---|---|---|---|---|---|---|
| 684.9821 | 560 | 630 | 740 | 209.2692 | 370 | 1911 | 7.673263 | 2.350456 | 30.55% |
df |> filter(Zone == "E") |>
summarise(Mean = mean(Income),
"1st Qu." = quantile(Income, probs = 0.25),
Median = median(Income),
"3st Qu." = quantile(Income, probs = 0.75),
Desv.pd = sd(Income),
Min = min(Income),
Max = max(Income),
Kurt = kurtosis(Income),
Skew = skewness(Income),
CV = sprintf("%.2f%%", (sd(Income)/mean(Income))*100)) |>
knitr::kable(caption = 'Medidas descritivas da variável Income (Renda)
das empresas da zona E')| Mean | 1st Qu. | Median | 3st Qu. | Desv.pd | Min | Max | Kurt | Skew | CV |
|---|---|---|---|---|---|---|---|---|---|
| 767.3879 | 594 | 705 | 864 | 242.41 | 470 | 1860 | 3.478582 | 1.689386 | 31.59% |
df |> ggplot(aes(Zone, Income)) +
stat_boxplot(geom = 'errorbar', na.rm = T, width = .25) +
geom_boxplot(na.rm = T) +
theme_test()
As zonas A e E são as que possuem a maior quantidade de empresas consideradas grandes, 30 e 23 respectivamente. Ou seja, juntas, essas duas zonas possuem 66,25% do total de empresas grandes. Consequentemente, elas possuem as maiores rendas médias, a zona A possui uma renda média de 652,2834 e a zona E possui uma renda média de 767,3879.
Por outro lado, as zonas que possuem as menores rendas médias são a zonas B e C, respectivamente. A zona B possui uma renda média de 320,7469 e a zona C possui uma renda média de 331,0195. O fato da zona C ter apenas uma empresa considerada grande, mas uma renda média maior do que a zona B que possui 13 é explicada pelo coeficiente de variação, a zona B possui uma variação maior da renda em relação à zona C.
Além disso, graficamente, é possível visualizar através do boxplot que todas as zonas possuem valores discrepantes, entretanto a zona B chama a atenção por possuir uma quantidade maior do que as outras zonas. Isso é evidenciado pelo seu coeficiente de variação de 70,86%, valor bem mais alto do que as demais.
1.5 Renda por SPAM
df |> ggplot(aes(SPAM, Income)) +
stat_boxplot(geom = 'errorbar', na.rm = T, width = .25) +
geom_boxplot(na.rm = T) +
theme_test()
# renda por empresas que se autopromoveram
df |> filter(SPAM == "yes") |>
summarise(Mean = mean(Income),
"1st Qu." = quantile(Income, probs = 0.25),
Median = median(Income),
"3st Qu." = quantile(Income, probs = 0.75),
Desv.pd = sd(Income),
Min = min(Income),
Max = max(Income),
Kurt = kurtosis(Income),
Skew = skewness(Income),
CV = sprintf("%.2f%%", (sd(Income)/mean(Income))*100)) |>
knitr::kable(caption = 'Medidas descritivas da variável Income (Renda)
das empresas que se automproveram na internet')| Mean | 1st Qu. | Median | 3st Qu. | Desv.pd | Min | Max | Kurt | Skew | CV |
|---|---|---|---|---|---|---|---|---|---|
| 436.782 | 233.5 | 389 | 578.5 | 275.8525 | 1 | 2510 | 4.390324 | 1.499564 | 63.16% |
df |> filter(SPAM == "no") |>
summarise(Mean = mean(Income),
"1st Qu." = quantile(Income, probs = 0.25),
Median = median(Income),
"3st Qu." = quantile(Income, probs = 0.75),
Desv.pd = sd(Income),
Min = min(Income),
Max = max(Income),
Kurt = kurtosis(Income),
Skew = skewness(Income),
CV = sprintf("%.2f%%", (sd(Income)/mean(Income))*100)) |>
knitr::kable(caption = 'Medidas descritivas da variável Income (Renda)
das empresas que não se autopromoveram na internet')| Mean | 1st Qu. | Median | 3st Qu. | Desv.pd | Min | Max | Kurt | Skew | CV |
|---|---|---|---|---|---|---|---|---|---|
| 424.7086 | 229 | 392 | 570 | 252.5118 | 1 | 1620 | 2.057913 | 1.120004 | 59.46% |
Não há grande diferença na renda média das empresas que se autopromoveram e das empresas que não se autopromoveram com propagandas na internet. Além disso, as outras medidas descritivas são bem próximas.
2 2. Amostragem Bernoulli
Sobre uma amostragem Bernoulli, temos que o estimador de Horvitz-Thompsom (\(\pi\)-estimador) para o total populacional é dado pela seguinte expressão:
\[ \hat{t}_\pi = \frac{1}{\pi}\sum_{k \in s}y_k, \] onde \(\pi, \pi \in [0,1]\) é a probabilidade de inclusão de primeira ordem, e \(y_k\) é observação que foi selecionada por meio do algoritmo de seleção. Além disso, temos que a variância estimada do (\(\pi\)-estimador) é dada por:
\[ \hat{V}(\hat{t}_\pi) = \frac{1}{\pi} \left( \frac{1}{\pi} - 1\right)\sum_{k \in s}y_k^{2}. \]
Agora, iremos fixar \(\pi=0,2\) e estimar o total da variável Income (Renda) e sua variância estimada. Primeiramente, utilizaremos a semente 13032024, em seguida, por meio da função S.BE() do pacote TeachingSampling, fixamos a probabilidade de inclusão de primeira ordem \(\pi = 0.2\) e selecionamos uma amostra a partir da nossa população. Após esses passos, obtivemos uma amostra de tamanho 455 do dataset Lucy. A função E.BE(), também do pacote TeachingSampling nos retorna na coluna y o valor estimado para o \(\pi\)-estimador, o erro padrão estimado e o coeficiente de variação estimado.
set.seed(13032024)
pi = 0.2
N = dim(df)[1] # tamanho da amostra do dataset Lucy
amostra = S.BE(N, pi) # amostra bernoulli selecionada
amostra = amostra[amostra != 0]
df_be = Lucy[amostra,]
E.BE((Lucy[amostra,]$Income), pi) |> knitr::kable(digits = 4)| N | y | |
|---|---|---|
| Estimation | 2275.0000 | 975810.0000 |
| Standard Error | 95.3939 | 48132.7772 |
| CVE | 4.1931 | 4.9326 |
| DEFF | Inf | 3.5976 |
Calculando manualmente, o estimador de Horvitz-Thompson (\(\pi\)-estimador) e sua respectiva variância estimada, temos os seguintes resultados:
# estimador de horvitz-thompsom para a variavel Income(Renda)
HT = (1/pi)*sum(df_be$Income);HT [1] 975810# variancia estimada do estimador de horvitz-thompsom
V_est_HT = (1/pi)*((1/pi)-1)*sum((df_be$Income)^2);V_est_HT[1] 23167642403 3. Tamanho Amostral Esperado
O tamanho amostral em uma amostragem Bernoulli é uma variável aleatória com distribuição Binomial de parâmetros \(N\) e \(\pi\). Portanto, o tamanho amostral esperado é:
\[ E[n(s)] = N\pi \]
N*pi[1] 479.2O tamanho amostral obtido no item anterior foi de \(n = 455\), já o tamanho amostral esperado foi de 479,2. Dado que o valor esperado leva em consideração todas as possíveis amostras, o que na prática é quase impossível, os valores do tamanho amostral esperado e do tamanho amostral observado são próximos.
4 4. Amostragem Bernoulli
Agora, iremos realizar uma amostragem Bernoulli considerando cada categoria da variável Level (Big, Medium e Small) como sendo uma população diferente. Para a categorias Big utilizaremos \(\pi = 0,2\), para a categoria Medium utilizaremos \(\pi = 0,30\) e para a categoria Small utilizaremos \(\pi = 0,25\). Em cada caso, será calculado o \(\pi\)-estimador e sua respectiva variância estimada.
4.1 Variável Level - Categoria Big
set.seed(13032024)
dfBig = df |> filter(Level == "Big")
pi_B = 0.20
NBig = dim(dfBig)[1]
amostra_B = S.BE(NBig, pi_B)
amostraBE = amostra_B[amostra_B != 0]
dfBBE = dfBig[amostraBE,]
E.BE((dfBBE$Income), pi_B) |> knitr::kable(digits = 4)| N | y | |
|---|---|---|
| Estimation | 60.0000 | 79560.0000 |
| Standard Error | 15.4919 | 20944.7760 |
| CVE | 25.8199 | 26.3258 |
| DEFF | Inf | 24.0832 |
Calculando manualmente, o estimador de Horvitz-Thompson (\(\pi\)-estimador) e sua respectiva variância estimada da categoria Big obtivemos os seguintes resultados:
# estimador de horvitz-thompsom para a variavel Income(Renda)
HT_Big = (1/pi_B)*sum(dfBBE$Income);HT_Big[1] 79560# variancia estimada do estimador de horvitz-thompsom
V_est_HT_Big = (1/pi_B)*((1/pi_B)-1)*sum((dfBBE$Income)^2);V_est_HT_Big[1] 4386836404.2 Variável Level - Categoria Medium
set.seed(13032024)
dfMedium = df |> filter(Level == "Medium")
pi_M = 0.30
NMedium = dim(dfMedium)[1]
amostra_M = S.BE(NMedium, pi_M)
amostraMBE = amostra_M[amostra_M != 0]
dfMBE = dfMedium[amostraMBE,]
E.BE((dfMBE$Income), pi_M) |> knitr::kable(digits = 4)| N | y | |
|---|---|---|
| Estimation | 713.3333 | 468633.3333 |
| Standard Error | 40.7976 | 27315.4495 |
| CVE | 5.7193 | 5.8287 |
| DEFF | Inf | 26.7527 |
Calculando manualmente, o estimador de Horvitz-Thompson (\(\pi\)-estimador) e sua respectiva variância estimada da categoria Medium obtivemos os seguintes resultados:
# estimador de horvitz-thompsom para a variavel Income(Renda)
HT_Medium = (1/pi_M)*sum(dfMBE$Income);HT_Medium[1] 468633.3# variancia estimada do estimador de horvitz-thompsom
V_est_HTMedium = (1/pi_M)*((1/pi_M)-1)*sum((dfMBE$Income)^2);V_est_HTMedium[1] 7461337804.3 Variável Level - Categoria Small
set.seed(13032024)
dfSmall = df |> filter(Level == "Small")
pi_S = 0.25
NSmall = dim(dfSmall)[1]
amostra_S = S.BE(NSmall, pi_S)
amostraSBE = amostra_S[amostra_S != 0]
dfSBE = dfSmall[amostraSBE,]
E.BE((dfSBE$Income), pi_S) |> knitr::kable(digits = 4)| N | y | |
|---|---|---|
| Estimation | 1544.0000 | 432668.0000 |
| Standard Error | 68.0588 | 20809.5060 |
| CVE | 4.4080 | 4.8096 |
| DEFF | Inf | 6.2324 |
Calculando manualmente, o estimador de Horvitz-Thompson (\(\pi\)-estimador) e sua respectiva variância estimada da categoria Small obtivemos os seguintes resultados:
# estimador de horvitz-thompsom para a variavel Income(Renda)
HT_Small = (1/pi_S)*sum(dfSBE$Income);HT_Small[1] 432668# variancia estimada do estimador de horvitz-thompsom
V_est_HTSmall = (1/pi_S)*((1/pi_S)-1)*sum((dfSBE$Income)^2);V_est_HTSmall[1] 4330355405 5. Sujestões para uma Estratégia Amostral
A partir da análise descritiva realizada, é possível perceber as variações de renda, tanto por zonas, quanto por dimensão da empresa. Nesse sentido, uma estratégia amostral eficiente deve levar em consideração esses aspectos. Portanto, uma sujestão seria um desenho amostral estratificado, visto que existem subcategorias dentro da população de interesse.
Assim, dependendo do obejtivo de estudo, os extratos poderiam ser as zonas (A, B, C, D e E) da variável Zone (Zona) ou as dimensões das empresas (Small, Medium e Big) da variável Level (Nível). Consequentemente, seria dado um tratamento particular a cada subgrupo da população, algumas vantagens seriam:
Representatividade da amostra selecionada em relação à população de interesse. Ao selecionar uma amostra aleatória simples da população, poderia acontecer de nenhuma empresa grande fosse selecionada, por exemplo.
Redução da variância na estimativa. Empresas do mesmo nível possuem rendas próximas, e empresas de nível diferentes possuem rendas distintas. A variância é reduzida porque os estratos são homogêneos internamente, mas heterogêneos entre si.