Теорија узорака

7. Из популације са вредностима обележја {1,2,…,100} извадити 15 простих случајних узорака без понављања обима 20. На сваком од њих наћи оцену укупне суме обележја и испитати који од узорака је најрепрезентативнији, тј. где је реализована вредност статистике \(T = N\overline{X_n}\) најближа стварној вредности.

vrednosti_obelezja = 1:100 
N = length(vrednosti_obelezja) # obim populacije
k = 15 # broj uzoraka
n = 20 # obimi uzoraka

# 15 puta vadimo psu bp obima n, smestamo ih u listu
uzorci = list()
for(i in 1:15) {
  uzorci[[i]] = sample(vrednosti_obelezja, n, replace = F)
}
uzorci

## [[1]]
##  [1] 92 57 18 47 62 60 30 43 58 41 71  1 44 54 83  5 93 19 39 32
## 
## [[2]]
##  [1]  48  99  41   7  39  70  59  65   8  58  57  45  27 100  22  38  18  74  28
## [20]  49
## 
## [[3]]
##  [1]  25  50  16  92  46  60  71  80  64  83  21  17  14  72  45  91  20  18  19
## [20] 100
## 
## [[4]]
##  [1] 93 53  6 72 83 95 46 25 77 98 67 48 13 36 30 31 73 85 92 34
## 
## [[5]]
##  [1] 98 94 69 85 72 33 82 10 16 24  3 23  2 61 59 27 73 62 68 18
## 
## [[6]]
##  [1]  73  54  93  32   5  85   8  76  89  96 100  68  60  15   1  80  46  31  61
## [20]  11
## 
## [[7]]
##  [1] 57 11 60 70 83 67 42 62 44 68 95 61 54 90 37 98 29 87  2 33
## 
## [[8]]
##  [1]  36  73  48  69  81  60  46  56  29  99 100   7  67  68  89  94  38  41  59
## [20]  20
## 
## [[9]]
##  [1] 83 18 74 20 13 39 64 33 26  1 24 89 87 51 63 34 50 46 49 55
## 
## [[10]]
##  [1] 39 57 54 60 50 15 95 96 38 90 66 30 51 93 87 45 13  7 11 53
## 
## [[11]]
##  [1] 53 45 74 38 29 30 12 72 10 83 76 35 52 49 69 28 63 20 36  2
## 
## [[12]]
##  [1] 17 43 47 62 76 86 79 31 49 61 73 55 64 12 26  7 29 84 38 83
## 
## [[13]]
##  [1] 58 56  7 98 45 91 85 72 67 71 93 92 79 76 46 52  8 36 26 16
## 
## [[14]]
##  [1]  61  44  65  67  97  95  54  16  48  47 100  30  71  31  26  58  73  94   2
## [20]  25
## 
## [[15]]
##  [1]  97  38  88  53  27  72   2   8  61  64  29  49  79  91  52  65  51  14 100
## [20]  55

# realizovane vrednosti statistike t_ocena = N*Xn na uzorcima
t_ocene = c()
for(i in 1:length(uzorci)) {
  t_ocene[i] = N*mean(uzorci[[i]])
}
t_ocene

##  [1] 4745 4760 5020 5785 4895 5420 5750 5900 4595 5250 4380 5110 5870 5520 5475

# suma obelezja na populaciji
t = sum(vrednosti_obelezja)
t

## [1] 5050

# which.min(abs(t_ocene - t)) - vraca indeks vrednosti iz vektora t_ocene gde je apsoulutna vrednost razlike t_ocene-t najmanja
t_najbliza = t_ocene[which.min(abs(t_ocene - t))]
t_najbliza

## [1] 5020

# which(t_ocene == t_najbliza) - vraca indeks vrednosti iz vektora t_ocene koja je jednaka vrednosti t_najbliza
najreprezentativniji = which(t_ocene == t_najbliza)
najreprezentativniji_uzorak = uzorci[[najreprezentativniji]]
najreprezentativniji_uzorak

##  [1]  25  50  16  92  46  60  71  80  64  83  21  17  14  72  45  91  20  18  19
## [20] 100

Прост случајан узорак

4. У кутији се налази 100 папирића на којима се налазе каро, пик, херц и треф. Десет пута се извлачи папирић, забележи се знак, а затим се папирић врати у кутију. Добијен је узорак (К, П, П, Х, К, Т, Т, П, Х, Х).

Ако каро има вредност 1 динар, пик 2, херц 3, а треф 4, оценити укупну суму новца која се налази у кутији.

б) Одредити оцену броја пикова.

в) Одредити непристрасне оцене дисперзија ових оцена.

N = 100
x = 1:4 # karo 1, pik 2, herc 3, tref 4

uzorak = c(1, 2, 2, 3, 1, 4, 4, 2, 3, 3)
n = length(uzorak)

t_ocena = N*mean(uzorak)
t_ocena

## [1] 250

# A - broj pikova od 100 papirica
br_pikova_u_uzorku = sum(uzorak==2)
pn = br_pikova_u_uzorku/n
A_ocena = N*pn
A_ocena

## [1] 30

D_t_ocena = N^2*var(uzorak)/n
D_t_ocena

## [1] 1166.667

D_A_ocena = N^2*pn*(1-pn)/(n-1)
D_A_ocena

## [1] 233.3333

5. За који од следећих планова простог слуучајног узорковања без понављања ће бити добијена најпрецизнија оцена средње вредности обележја на популацији? Претпоставимо да обележје на свакој од популација има дисперзију 100.

а) Узорак обима 400, добијен из популације обима 4000.

б) Узорак обима 30, добијен из популације обима 200.

в) Узорак обима 3000, добијен из популације обима 300000000.

# trazimo gde je D(Xn) najmanja

s2 = 100

D_xn = function(N,n) {
  (N-n)*s2/(N*n)
}

D_xn_1 = D_xn(4000, 400)
D_xn_2 = D_xn(200, 30)
D_xn_3 = D_xn(300000000, 3000)
D = c(D_xn_1, D_xn_2, D_xn_3)
which.min(D)

## [1] 3

8. Оцењује се број дивљих животиња у неком региону који је подељен на 286 области. Изабран је прост случајан узорак без понављања од 15 области и дат је број животиња у њима: 1, 50, 21, 98, 2, 36, 4, 29, 7, 15, 86, 10, 21, 5, 4.

а) Наћи узорачку дисперзију, оцену просечног броја животиња по области, оцену дисперзије те оцене, као и оцену стандардне грешке.

б) Наћи 90% интервал поверења за укупан број животиња.

в) Коју величину узорка треба узети за оцењивање укупног броја животиња тако да одступање не буде веће од 2000 животиња са вероватноћом 0.9.

N = 286
n = 15
uzorak = c(1, 50, 21, 98, 2, 36, 4, 29, 7, 15, 86, 10, 21, 5, 4)

sn2 = var(uzorak)
sn2

## [1] 919.0667

xn = mean(uzorak)
xn

## [1] 25.93333

D_xn_ocena = (N-n)*sn2/(N*n)
D_xn_ocena

## [1] 58.05759

# xn je nepristrasna ocena, pa je MSE(xn) = D(xn), pa je standardna greska
# koren iz disperzije, i ocena standardne greske = koren iz ocene disperzije
st_greska_xn_ocena = sqrt(D_xn_ocena)
st_greska_xn_ocena

## [1] 7.619553

# n = 15 < 30 ---> Studentov kvantil
alpha = 0.1
z_student = qt(1 - alpha/2, n-1)
z_student

## [1] 1.76131

t_ocena = N*mean(uzorak)
I_poverenja_za_t_90 = c(t_ocena - z_student*sqrt((N-n)*N*var(uzorak)/n),
                        t_ocena + z_student*sqrt((N-n)*N*var(uzorak)/n))
I_poverenja_za_t_90

## [1]  3578.70 11255.17

d = 2000
alpha = 0.1
z = qnorm(1 - alpha/2)
n1 = (d^2/(z^2*sn2*N^2) + 1/N)^{-1}
n = ceiling(n1)
n

## [1] 44

10. Посматрамо популацију обима 15 са вредностима обележја 3, 4, 12, 9, 8, 2, 6, 15, 14, 17, 1, 7, 12, 0, 1 и 5. Размотримо принцип простог случајног узорковања без понављања за узорак обима 4.

а) Одредити све могуће узорке и вероватноће да сваки од тих узорака буде изабран.

б) За сваки узорак израчунати узорачку средњу вредност, узорачку дисперзију, оцену дисперзије узорачке средине и медијану.

в) Показати да је средња вредност обележја на узорку непристрасна оцена средње вредности обележја популације и да је узорачка дисперзија непристрасна оцена дисперзије обележја на популацији.

obelezje = c(3, 4, 12, 9, 8, 2, 6, 15, 14, 17, 1, 7, 12)
N = length(obelezje)
n = 4

br_uzoraka = choose(N, n)
br_uzoraka

## [1] 715

uzorak = list()
s = 1
for(i in 1:N) {
  for(j in 1:N) {
    for(k in 1:N) {
      for(l in 1:N) {
      if((i < j) & (j < k) & (k < l)) {
        uzorak[[s]] = c(obelezje[i], obelezje[j],
                        obelezje[k], obelezje[l])
        s = s+1
        }
      }
    }
  }
}
length(uzorak)

## [1] 715

# verovatnoca svakog od uzoraka je jednaka 1/br_uzoraka

xn = c()
sn2 = c()
D_xn_ocena = c()
medijana = c()
for(i in 1:br_uzoraka) {
  xn[i] = mean(uzorak[[i]])
  sn2[i] = var(uzorak[[i]])
  D_xn_ocena[i] = (N-n)*sn2[i]/(N*n)
  medijana[i] = median(uzorak[[i]])
}
xn[1:10] # prikazujemo prvih 10

##  [1] 7.00 6.75 5.25 6.25 8.50 8.25 9.00 5.00 6.50 7.75

sn2[1:10]

##  [1] 18.00000 16.91667 20.91667 16.25000 35.00000 30.91667 44.66667 23.33333
##  [9] 16.33333 24.25000

D_xn_ocena[1:10]

##  [1] 3.115385 2.927885 3.620192 2.812500 6.057692 5.350962 7.730769 4.038462
##  [9] 2.826923 4.197115

medijana[1:10]

##  [1] 6.5 6.0 3.5 5.0 8.0 8.0 8.0 3.5 5.5 8.0

# nepristrasnost
Exn = sum(xn * 1/br_uzoraka)
Exn

## [1] 8.461538

mean(obelezje)

## [1] 8.461538

Esn2 = sum(sn2 * 1/br_uzoraka)
Esn2

## [1] 27.26923

var(obelezje)

## [1] 27.26923

Узорковање са неједнаким вероватноћама

3. Популација се састоји од 15 јединки чије су величине 23, 30, 41, 26, 53, 60, 28, 52, 113, 72, 80, 35, 42, 38 и 52. Изабрати узорак са понављањем и узорак без понављања обима 5 са вероватноћама пропорционалним величинама користећи Лахиријев метод.

velicine = c(32, 21, 24, 81, 29, 22, 62, 3, 56, 65, 100, 87, 68, 83, 77, 51, 76, 67, 38, 97)
N = length(velicine)
sum(velicine)

## [1] 1139

n = 15
K = max(velicine)
K

## [1] 100

# u svakom koraku biramo uredjeni par (i, R)
# ako R <= Mi - biramo i-tu jedinku u uzorak
# ako R > Mi - nista, nastavljamo postupak
# uzorak sp - dok ne izaberemo n
# uzorak bp - dok ne izaberemo n razlicitih
uzorak_sp = c()
while(length(uzorak_sp) < n) {
  i = sample(1:N, 1) 
  R = sample(1:K, 1)
  if(R <= velicine[i]) {
    uzorak_sp = c(uzorak_sp, i)
  }
}
uzorak_sp

##  [1] 11 20 11 12 11 18 20  9 18 10  9  1 12 12 13

uzorak_bp = c()
while(length(uzorak_bp) < n) {
  i = sample(1:N, 1)
  R = sample(1:K, 1)
  if(R <= velicine[i] & !(i %in% uzorak_bp)) {
    uzorak_bp = c(uzorak_bp, i)
  }
}
uzorak_bp

##  [1]  4 15 16  5 11 18 12 10 13 17  1 20 19  3 14

7. У датотеци populacija налазе се величине 120 јединки из популације. Извадити узорак са понављањем и узорак без понављања обима 90 користећи метод кумуланти.

populacija = read.table("populacija.txt")
M = sum(populacija$velicine)
n = 90

# formiramo kumulante
kumulante = cumsum(velicine)
kumulante = c(0, kumulante)

# biramo slucajan broj B: 0 < B <= M
# ako M1 + ... + Mi < B <= M1 + ... + M{i+1} biramo jedinku (i+1) u uzorak
# uzorak sa ponavljanjem - sprovodimo postupak n puta
uzorak = c()
while(length(uzorak) < n) {
  B = sample(0:M, 1)
  for(i in 1:(length(kumulante) - 1)) {
    if((kumulante[i] < B) & (B <= kumulante[i+1])) {
      uzorak = c(uzorak, i)
      break
    }
  }
}
uzorak

##  [1]  21 105  77  22 116  67  21   6  98  38  59  67  22  16  82  20  84  77  58
## [20]  33 109  16 114  70  35 100  24  79  60  88  71 104  81   9  66  16 116  90
## [39]  23  26 105  58  34  20   6 119 103  79  34  18  98  38   1  25  18  23  63
## [58]   1  29 112  21  12  82  92 120  53  36  59  64 102  93  15  67  46  95  86
## [77]  53  26 102  44 106 117  51  20   1  82  99  53  34  63

length(unique(uzorak)) # broj razlicitih jedinki u uzorku

## [1] 59

# uzorak bez ponavljanja - sprovodimo postupak dok ne dobijemo n razlicitih
uzorak = c()
while (length(uzorak) < n) {
  B = sample(0:M, 1)
  for(i in 1:(length(kumulante) - 1)) {
    if((kumulante[i]<B) & (B<=kumulante[i+1]) & (!(i %in% uzorak))) {
      uzorak = c(uzorak, i)
      break
    }
  }
}
uzorak

##  [1]  13 100 114  51  95  85 112  33  53 118  74  69  21  32   8 119  54  72  73
## [20] 113  23   6  16 101  10  39  18 120  48  98  36  77 105  20  56  67  61   9
## [39]  75  94  12  84  63  58  14 109  15  71  66   3  17 108  46  19  70 106  62
## [58]  93 115  78  89 104  26  29 111  35  38  99  22  59  81  31  37  11   1  87
## [77] 117  30  40   7  43  47  82  90  60 102  86  64  88  91

length(unique(uzorak))

## [1] 90

8. Из популације коју чине три поља на којима се узгаја пшеница бира се узорак обима 2 са вероватноћама пропорционалним величинама, са понављањем. Дати су подаци о количини произведене пшенице на сваком пољу и вероватноће избора сваког поља.

i: 1 | 2 | 3
x_i: 11 | 6 | 25
p_i: 0.3 | 0.2 | 0.5

Налажењем свих могућих узорака одредити вероватноће укључења сваког елемента, као и Хансен-Хурвицову и Хорвиц-Томпсонову оцену за укупну производњу пшенице за сваки узорак.

N = 3
n = 2

i = 1:3
xi = c(11, 6, 25)
pi = c(0.3, 0.2, 0.5)

uzorci = list()
obelezje_na_uzorku = list()
i = 1
for (j in 1:N) {
  for (k in 1:N) {
    uzorci[[i]] = c(j, k)
    obelezje_na_uzorku[[i]] = c(xi[j], xi[k])
    i = i+1
  }
}
length(uzorci)

## [1] 9

N^n

## [1] 9

# v-ca uzorka (1,1) je p1*p1, (1,2) p1*p2, itd.
p_uzorka = c()
t_hh = c()
for(i in 1:length(uzorci)) {
  p_uzorka[i] = prod(pi[uzorci[[i]]])
  t_hh[i] = sum(xi[uzorci[[i]]]/pi[uzorci[[i]]])/n
}
p_uzorka

## [1] 0.09 0.06 0.15 0.06 0.04 0.10 0.15 0.10 0.25

sum(p_uzorka)

## [1] 1

# realizovane vrednosti ocene t_hh na svim mogucim uzorcima
t_hh

## [1] 36.66667 33.33333 43.33333 33.33333 30.00000 40.00000 43.33333 40.00000
## [9] 50.00000

# verovatnoce ukljucenja
pii = 1 - (1-pi)^n
pii

## [1] 0.51 0.36 0.75

# ili zbir verovatnoca svih uzoraka koji sadrze i 
pii1 = rep(0, N)
for(i in 1:N) { 
  for(s in 1:length(uzorci)) {
    if(i %in% uzorci[[s]]) {
      pii1[i] = pii1[i] + p_uzorka[s]
    }
  }
}
pii1

## [1] 0.51 0.36 0.75

pii

## [1] 0.51 0.36 0.75

# realizovane vrednosti ocene t_ht
t_ht = c()
for(i in 1:length(uzorci)) {
  if(uzorci[[i]][1] == uzorci[[i]][2]) { # ako ista jedinka odabrana 2 puta
    t_ht[i] = obelezje_na_uzorku[[i]][1]/pii[uzorci[[i]][1]]
  }
  else { # ako dve razlicite jedinke u uzorku
    t_ht[i] = sum(obelezje_na_uzorku[[i]]/pii[uzorci[[i]]])
  }
}
t_ht

## [1] 21.56863 38.23529 54.90196 38.23529 16.66667 50.00000 54.90196 50.00000
## [9] 33.33333

# provera - ispitujemo nepristrasnost
E_t_hh = sum(t_hh*p_uzorka)
E_t_hh

## [1] 42

E_t_ht = sum(t_ht*p_uzorka)
E_t_ht

## [1] 42

sum(xi)

## [1] 42

10. Дата је популација од четири прашуме, њихове површине и бројеви тигрова који живе у њима:

Прашума: 1 | 2 | 3 | 4

Површина (у km2): 100 | 200 | 300 | 500

Број тигрова: 11 | 20 | 23 | 54

Вади се узорак од 2 прашуме. Оценити укупан број тигрова, одредити дисперзију те оцене и наћи непристрасну оцену те дисперзије ако је узорковање вршено:

а) Са понављањем (користити Хансен-Хурвицову оцену),

б) Без понављања (користити Хорвиц-Томпсонову оцену),

и извучен је узорак (1,2).

N = 4
n = 2

br_tigrova = c(11, 20, 23, 54)
Mi = c(100, 200, 300, 500)
M = sum(Mi)

uzorak = br_tigrova[c(1, 2)]

pi = Mi/M
pi

## [1] 0.09090909 0.18181818 0.27272727 0.45454545

sum(pi)

## [1] 1

# sa ponavljanjem
t_hh = sum(uzorak/pi[c(1, 2)])/n
t_hh

## [1] 115.5

t = sum(br_tigrova)
t

## [1] 108

D_t_hh = sum(pi*(br_tigrova/pi - t)^2)/n
D_t_hh

## [1] 110.9333

D_t_hh_ocena = sum((uzorak/pi[c(1, 2)] - t_hh)^2)/(n*(n-1))
D_t_hh_ocena

## [1] 30.25

# bez ponavljanja
# za HT ocenu potrebne su nam verovatnoce ukljucenja
# uzorak bp -> ne vazi pii = (1 - pi)^n u opstem slucaju!
# pii = suma verovatnoca svih uzoraka koji sadrze i-tu jedinku
uzorci = list()
i = 1
for(j in 1:N) {
  for(k in 1:N) {
    if(j < k) {
      uzorci[[i]] = c(j, k)
      i = i + 1
    }
  }
}

# v-ca uzorka (i,j) =
# P(odaberemo i)P(odaberemo j|odabrali smo i) + P(odaberemo j)P(odaberemo i|odabrali smo j) =
# (Mi/M)*(Mj/(M-Mi)) + (Mj/M)*(Mi/(M-Mj)) 
vca_uzorka = function(i, j) {
  (Mi[i]/M)*(Mi[j]/(M - Mi[i])) + (Mi[j]/M)*(Mi[i]/(M - Mi[j]))
}

p_uzorka = c()
for(i in 1:length(uzorci)) {
  p_uzorka[i] = vca_uzorka(uzorci[[i]][1], uzorci[[i]][2])
}
p_uzorka

## [1] 0.03838384 0.06136364 0.12121212 0.12878788 0.25252525 0.39772727

sum(p_uzorka)

## [1] 1

pii = rep(0,N)
for(i in 1:N) {
  for(j in 1:length(uzorci)) {
    if(i %in% uzorci[[j]]) {
      pii[i] = pii[i] + p_uzorka[j]
    }
  }
}
pii

## [1] 0.2209596 0.4196970 0.5878788 0.7714646

# provera - sum(pii) == n, za uzorak bez ponavljanja
# (nekad moze da vrati FALSE iako zapravo jeste =, zbog zaokruzivanja)
sum(pii) == n

## [1] TRUE

t_ht = sum(uzorak/pii[c(1, 2)])
t_ht

## [1] 97.43629

sum(br_tigrova)

## [1] 108

D_t_ht = sum((1 - pii)*br_tigrova^2/pii)
for(i in 1:N) {
  for(j in 1:N) {
    if(i != j) {
      pij = vca_uzorka(i, j)
      D_t_ht = D_t_ht + (pij - pii[i]*pii[j])*(br_tigrova[i]*br_tigrova[j])/(pii[i]*pii[j])
    }
  }
}
D_t_ht

## [1] 125.4964

D_t_ht_ocena = sum((1 - pii[c(1, 2)])*uzorak^2/(pii[c(1, 2)]^2))
for(i in c(1, 2)) {
  for(j in c(1, 2)) {
    if(i != j) {
      pij = vca_uzorka(i, j)
      D_t_ht_ocena = D_t_ht_ocena + (pij - pii[i]*pii[j])*(uzorak[i]*uzorak[j])/(pii[i]*pii[j]*pij)
    }
  }
}
D_t_ht_ocena

## [1] -3470.008

D_t_ht_ocena_syg = 0
for (i in c(1, 2)) {
  for(j in c(1, 2)) {
    if(i < j) {
      pij = vca_uzorka(i, j)
      D_t_ht_ocena_syg = D_t_ht_ocena_syg + (pii[i]*pii[j] - pij)*((uzorak[i]/pii[i] - uzorak[j]/pii[j])^2)/pij
    }
  }
}
D_t_ht_ocena_syg

## [1] 6.420883

Стратификовано узорковање

1. Испитује се колико су метара претрчали ученици једног разреда на часу физичког васпитања. У разреду има 112 ученика, од чега 59 девојчица и 53 дечака. Посматрано је 7 дечака и 8 девојчица. Забележено је:

Дечаци: 879, 810, 789, 567, 900, 870, 777.

Девојчице: 450, 234, 679, 456, 239, 555, 560, 467.

Оценити колико су просечно метара претрчала деца тог разреда и оценити дисперзију те оцене. Затим урадити исто коришћењем оцена предложених код простог случајног узорковања без `понављања.

N = 112
N_dec = 53
N_dev = 59

n_dec = 7
n_dev = 8

# broj pretrcanih metara za decake i devojcice
uzorak_decaci = c(879, 810, 789, 567, 900, 870, 777)
uzorak_devojcice = c(450, 234, 679, 456, 239, 555, 560, 467)

x_str = (N_dec*mean(uzorak_decaci) + N_dev*mean(uzorak_devojcice))/N
x_str

## [1] 617.7181

# ne mozemo da odredimo disperziju ocene jer nemamo vrednosti obelezja za celu populaciju
# racunamo nepristrasnu ocenu disperzije
D_x_str_ocena = (N_dec^2*var(uzorak_decaci)*(1-n_dec/N_dec)/n_dec 
                 + N_dev^2*var(uzorak_devojcice)*(1-n_dev/N_dev)/n_dev)/(N^2)
D_x_str_ocena

## [1] 1067.26

# PSU BP
n = n_dec + n_dev
uzorak = c(uzorak_decaci, uzorak_devojcice)
xsr_ocena = mean(uzorak)
xsr_ocena

## [1] 615.4667

D_xsr_ocena = var(uzorak)*(N-n)/(N*n)
D_xsr_ocena

## [1] 2822.814

D_x_str_ocena

## [1] 1067.26

3. У датотеци stratumi налазе се подаци о вредностима обележја на популацији која се састоји из 2 стратума (стратуми су различитих величина - да би могли да се сместе у две колоне базе, мањи стратум је допуњен NA вредностима до величине већег стратума. Избацити их при рачунању.). Бира се стратификован узорак обима 500, тако што се из мањег стратума бира прост случајан узорак са понављањем обима 200, а из већег обима 300.

а) Оценити укупну суму обележја популације, израчунати дисперзију те оцене и оцену те дисперзије.

stratumi = read.table("stratumi.txt")

# koliko ima NA vrednosti
sum(is.na(stratumi$str1))

## [1] 198

sum(is.na(stratumi$str2))

## [1] 0

strat = list()
# uzimamo samo vrednosti koje nisu NA
strat[[1]] = na.exclude(stratumi$str1) 
strat[[2]] = stratumi$str2

# velicine stratuma
Ni = sapply(strat,length)
Ni

## [1] 690 888

# obim populacije
N = sum(Ni)

n = 500
uzorak = list()
# vadimo PSU SP obima 200 iz prvog, a 300 iz drugog stratuma
uzorak[[1]] = sample(strat[[1]], 200, replace = T)
uzorak[[2]] = sample(strat[[2]], 300, replace = T)
uzorak

## [[1]]
##   [1] 477 933  48 484   1  68 795 103 945 735 376 749 125  33 976 209 600 376
##  [19] 368 246  64   9 895 945 651 465 763 275 708  23 274 413 524 438 895 125
##  [37] 719 338 747 341 646 794 868 452 902 801 980 226  48 423 933 604 367 338
##  [55] 737 439 793 101 536 234 696 654 227 394 649 934 441 653 307 263 241 422
##  [73] 885  15 380 632 404  64 474  10 860 397  36  64 544  74 438 477 238 846
##  [91] 529 274 113 604 270 674 980 377 895 491 332   9 438 729  29 202 508 112
## [109] 558  69 671 472 124 698 154 302 334 772 938 962 712  28 911 397 164 931
## [127] 483 719 657 934 945 460 312 318 114 334 536 205 513 457 957   6 860 651
## [145] 312 677 147 441 940  32 740 190 456 649 311 592 957 774  63  27 768 354
## [163] 227 397 792 965 151 367 449 349 226   1  89 700 124 679  35 768 517 530
## [181] 787 363 712 477 870 155 969 602 226 602 387 938 962 832 486 337 602 949
## [199] 956 341
## 
## [[2]]
##   [1] 366 196 865 923 536 477  39 590 513 401 621 598  48 740 577 385 714 535
##  [19] 699 387 756 931 869 713 877 577 411  19 401 786 360 594 203 474  48 366
##  [37] 827 108 156 746 656 776 411 865 170 101 668 411  40 647 289 903  54 236
##  [55] 668 292 417 674 350 499 993 358 496 636 310 146 278 339 496 875 664 240
##  [73]  99 898 228 835 715 168 992 608 934 541 169 549 886 916 951  58 394 940
##  [91] 209 594 974 257 951  46 320 639 302 365 547 636 438 437 448 376 667 715
## [109] 453 390  83 814 815 152 382 345 891 866 470 673 201 312  82 379  46 814
## [127] 699 511 368 109 107  46 609 658 790 617 652 537 498 271 879 882 119 107
## [145] 949 743 809 365 273 219 402 436 688 244 978 843 465 654 205 409 217 281
## [163]  72 668  28 651  11 663 219 281 241 134 696 162  72 271 315 430 324 198
## [181] 938 473 226 711 733  97 458 146 821 909 903 175 852 993 208 245 869 825
## [199] 581 577 289 499 314 374 125 475 950 639 105 778 143 416 201 353  14 511
## [217] 734 389 984 363 608 663 935 879 496 365 553  60 156 859 156 292 884 658
## [235] 170 538 458 836 181 618 487  79  44 430 108  51 859 569 278 110 147 289
## [253] 461 751 553 306 114  40 491 210 368 991 651 539 740 143 774 261 569 635
## [271] 927 203 711  82 156 416 830 594 395 780 611 851  46 905 733 362 360 756
## [289] 473 158 508 119 913 293 747 406 577 618 173 788

# obimi uzoraka
ni = sapply(uzorak, length)

# racunamo uzoracke srednje vrednosti
xni = sapply(uzorak, mean)
t_str = sum(Ni*xni)
t_str

## [1] 762282.4

# prava vrednost t
obelezje_na_populaciji = c(strat[[1]], strat[[2]])
sum(obelezje_na_populaciji)

## [1] 785903

# disperzija ocene t_str
s2i = sapply(strat, var)
D_t_str = sum(Ni*(Ni-1)*s2i/ni)
D_t_str

## [1] 420368534

# ocena disperzije t_str
sn2i = sapply(uzorak, var)
D_t_str_ocena = sum(Ni^2*sn2i/ni)
D_t_str_ocena

## [1] 412517409

4. а) Популација је подељена на три стратума чије су величине 123, 102 и 180, a одговарајуће дисперзије обележја су 116, 143 и 170, редом. Ако се бира стратификован случајан узорак обима 100, одредити величину узорка који се бира из сваког стратума користећи пропорционални избор.

б) Извести формулу за одређивање обима узорака који се узимају из сваког стратума ако се користи Нејманова метода, а затим одредити величину узорка који се узима из сваког стратума за популацију и дати обим узорка из дела а).

Ni = c(123, 102, 180)
si2 = c(116, 143, 170)
n = 100
N = sum(Ni)

# proporcionalni izbor
n_prop = round(n*Ni/N)
n_prop

## [1] 30 25 44

sum(n_prop) # mora biti = n (ako nije, neki smanjimo/povecamo tako da u zbiru daju n)

## [1] 99

while(sum(n_prop) > n) { # ako je > n, smanjujemo slucajno odabrani za 1 dok ne bude n
  pom = sample(1:length(n_prop), 1)
  n_prop[pom] = n_prop[pom] - 1
}
while(sum(n_prop) < n) { # ako je < n, povecavamo slucajno odabrani za 1 dok ne bude n
  pom = sample(1:length(n_prop), 1)
  n_prop[pom] = n_prop[pom] + 1
}
sum(n_prop)

## [1] 100

# Nejmanov izbor
n_nejman = Ni*sqrt(si2)*n / sum(Ni*sqrt(si2))
n_nejman = round(n_nejman)
sum(n_nejman)

## [1] 100

8. У датотеци netflix_titles налазе се подаци о филмовима и серијама додатим на платформу Netflix у периоду од 2008. до 2021. године.

Поделити наслове у стратуме по годинама када су додати, али тако да се сви наслови додати пре 2015. године налазе у једном стратуму. Изабрати стратификован случајан узорак без понављања обима 500 користећи Нејманов избор, па оценити удео серија у насловима.

б) Поделити серије у стратуме по рејтингу (постоји само једна серија са рејтингом TV-Y7-FV, њу сместити у стратум са серијама са рејтингом R). Изабрати стратификован случајан узорак без понављања обима 200 серија користећи Нејманов избор, па оценити просечан број сезона.

netflix_titles = read.csv("netflix_titles.csv")

# a)
sum(netflix_titles$date_added == "")

## [1] 10

sum(netflix_titles$type == "")

## [1] 0

indeksi = which(netflix_titles$date_added == "")
netflix_titles_full = netflix_titles[-indeksi, ] # izbacujemo opservacije gde nemamo datum dodavanja na netflix
sum(netflix_titles_full$date_added == "")

## [1] 0

N = length(netflix_titles_full$show_id)

netflix_titles_full$date_added[1]

## [1] "September 25, 2021"

# zelimo godinu -> koristimo strsplit da razdvojimo string
strsplit(netflix_titles_full$date_added[1], ", ")

## [[1]]
## [1] "September 25" "2021"

year_added = sapply(strsplit(netflix_titles_full$date_added, ", "), function(x) {x[2]}) # uzimamo godine
year_added = as.numeric(year_added)

netflix_titles_full = cbind(netflix_titles_full, year_added)

unique(netflix_titles_full$year_added)

##  [1] 2021 2020 2019 2018 2017 2016 2015 2014 2013 2012 2011 2009 2008 2010

# formiramo stratume
stratumi = list()
s = 1
for (i in unique(netflix_titles_full$year_added)[unique(netflix_titles_full$year_added) > 2014]) {
  stratumi[[s]] = netflix_titles_full[netflix_titles_full$year_added == i, ]
  s = s + 1
}
stratumi[[s]] = netflix_titles_full[netflix_titles_full$year_added <= 2014, ]

Ni = sapply(stratumi, function(x) {length(x[, 1])})
sum(Ni) == length(netflix_titles_full$show_id)

## [1] TRUE

n = 500
si = sapply(stratumi, function(x) {sd(x$type == "TV Show")})
ni = round(Ni*si*n / sum(Ni*si))
sum(ni)

## [1] 501

pom = sample(1:length(ni), 1)
ni[pom] = ni[pom] - 1
sum(ni)

## [1] 500

uzorak = list(); xni = c()
s = 1
for (i in 1:length(stratumi)) {
  uzorak[[s]] = stratumi[[s]][sample(1:Ni[s], ni[s]), ]
  xni[s] = mean(uzorak[[s]]$type == "TV Show")
}
(pstr = sum(Ni*xni) / N)

## [1] 0.3863636

# b)
TVS = netflix_titles[netflix_titles$type == "TV Show", ] # biramo samo TV serije

sum(TVS$rating == "")

## [1] 2

sum(TVS$duration == "")

## [1] 0

indeksi = which(TVS$rating == "")
TVS_full = TVS[-indeksi, ]
sum(TVS_full$rating == "")

## [1] 0

Ns = length(TVS_full$show_id)

# zelimo da izvucemo broj sezona iz Duration -> koristimo strsplit
TVS_full$duration[1]

## [1] "2 Seasons"

strsplit(TVS_full$duration[1], " ")

## [[1]]
## [1] "2"       "Seasons"

seasons = sapply(strsplit(TVS_full$duration, " "), function(x) {x[1]})
seasons = as.numeric(seasons)

TVS_full = cbind(TVS_full, seasons)
unique(TVS_full$rating)

## [1] "TV-MA"    "TV-14"    "TV-Y7"    "TV-PG"    "TV-Y"     "TV-G"     "R"       
## [8] "NR"       "TV-Y7-FV"

sum(TVS_full$rating == "TV-Y7-FV") # ne moze u zaseban stratum, samo 1 element

## [1] 1

# stavljamo u stratum sa "R"

stratumi = list()
s = 1
for (i in unique(TVS_full$rating)[unique(TVS_full$rating) != "TV-Y7-FV"]) {
  stratumi[[s]] = TVS_full[TVS_full$rating == i, ]
  s = s + 1
}

# ubacujemo opservacije sa rejtingom "TV-Y7-FV" u stratum sa rejtingom "R"
for (i in 1:length(stratumi)) {
  if (stratumi[[i]]$rating[1] == "R") { # svi imaju isti rejting, dovoljno je da pogledamo prvi element
    stratumi[[i]] = rbind(stratumi[[i]], TVS_full[TVS_full$rating == "TV-Y7-FV", ])
  }
}
Ni = sapply(stratumi, function(x) {length(x[, 1])})
sum(Ni) == Ns

## [1] TRUE

ns = 200
si = sapply(stratumi, function(x) {sd(x$seasons)})
ni = round(Ni*si*n / sum(Ni*si))
sum(ni)

## [1] 500

uzorak = list(); xni = c()
s = 1
for (i in 1:length(stratumi)) {
  uzorak[[s]] = stratumi[[s]][sample(1:Ni[s], ni[s]), ]
  xni[s] = mean(uzorak[[s]]$seasons)
}
(xstr = sum(Ni*xni) / N)

## [1] 0.5008372

Количничко оцењивање

3. Дата је популација обима 4 и вредности обележја које се испитује (x) и помоћног обележја (y):

x: 3, 5, 7, 9

y: 1, 2, 2, 3

На основу простих случајних узорака без понављања обима 2 и 3 испитати да ли је количничка оцена количника популације непристрасна.

N = 4
x = c(3, 5, 7, 9)
y = c(1, 2, 2, 3)

# formiramo sve uzorke bez ponavljanja obima 2
uzorci_2 = list()
i = 1 
for(j in 1:N) {
  for(k in 1:N) {
    if(j < k) {
      uzorci_2[[i]] = rbind(c(x[j], x[k]),
                            c(y[j], y[k])
                            )
      i = i + 1
    }
  }
}
uzorci_2

## [[1]]
##      [,1] [,2]
## [1,]    3    5
## [2,]    1    2
## 
## [[2]]
##      [,1] [,2]
## [1,]    3    7
## [2,]    1    2
## 
## [[3]]
##      [,1] [,2]
## [1,]    3    9
## [2,]    1    3
## 
## [[4]]
##      [,1] [,2]
## [1,]    5    7
## [2,]    2    2
## 
## [[5]]
##      [,1] [,2]
## [1,]    5    9
## [2,]    2    3
## 
## [[6]]
##      [,1] [,2]
## [1,]    7    9
## [2,]    2    3

# racunamo vrednosti koje uzima slucajna velicina R_ocena = xn/yn na uzorcima
# posto je u pitanju prosto slucajno uzorkovanje, svaki od uzoraka ima jednaku verovatnocu 1/length(uzorci_2)
# (1/broj uzoraka)
# Dakle, sl.velicina R_ocena uzima vrednosti iz vektora R_ocena_2, svaku sa istom verovatnocom 1/length(uzorci_2)
R_ocena_2 = c()
for(i in 1:length(uzorci_2)) {
  R_ocena_2[i] = mean(uzorci_2[[i]][1, ]/uzorci_2[[i]][2, ])
}

ER_ocena_2 = sum(R_ocena_2 * 1/length(uzorci_2))
ER_ocena_2

## [1] 3

# kolicnik obelezja populacije
R = mean(x)/mean(y)
R

## [1] 3

# ispostavlja se da jeste nepristrasna u ovom slucaju

# n = 3
uzorci_3 = list()
i = 1
for(j in 1:N) {
  for(k in 1:N) {
    for(s in 1:N) {
      if(j < k && k < s) {
        uzorci_3[[i]] = rbind(c(x[j],x[k],x[s]),
                              c(y[j],y[k],y[s])
                              )
        i = i + 1
      }
    }
  }
}
uzorci_3

## [[1]]
##      [,1] [,2] [,3]
## [1,]    3    5    7
## [2,]    1    2    2
## 
## [[2]]
##      [,1] [,2] [,3]
## [1,]    3    5    9
## [2,]    1    2    3
## 
## [[3]]
##      [,1] [,2] [,3]
## [1,]    3    7    9
## [2,]    1    2    3
## 
## [[4]]
##      [,1] [,2] [,3]
## [1,]    5    7    9
## [2,]    2    2    3

R_ocena_3 = c()
for(i in 1:length(uzorci_3)) {
  R_ocena_3[i] = mean(uzorci_3[[i]][1, ]/uzorci_3[[i]][2, ])
}

ER_ocena_3 = sum(R_ocena_3 * 1/length(uzorci_3))
ER_ocena_3

## [1] 3

# jeste nepristrasna

4. Дате су следеће вредности:

Школа | 1 | 2 | 3 | 4

Број деце која имају 5 из математике | 200 | 300 | 100 | 250

Број посматраних разреда | 4 | 5 | 3 | 4

Истраживач жели да посматрањем прве и треће школе донесе закључак о укупном броју деце која имају пет из математике. Испитати да ли је количничка оцена која се користи непристрасна (ако се користи просто случајно узорковање без понављања), одредити њену вредност на датом узорку, а затим упоредити са оценом добијеном на основу простог случајног узорка без понављања. Која је боља?

br_petica = c(200, 300, 100, 250)
br_razreda = c(4, 5, 3, 4)

# uzorak koji je odabran
indeks = c(1, 3)
uzorak_br_petica = br_petica[indeks]
uzorak_br_razreda = br_razreda[indeks]

# kolicnicka ocena na osnovu odabranog uzorka
R_ocena = mean(uzorak_br_petica) / mean(uzorak_br_razreda)
t_r = R_ocena * sum(br_razreda)
t_r

## [1] 685.7143

N = 4
# nepristrasnost kolicnicke ocene
uzorci = list()
i = 1
for(j in 1:N) {
  for(k in 1:N) {
    if(j < k) {
      uzorci[[i]] = rbind(c(br_petica[j], br_petica[k]),
                          c(br_razreda[j], br_razreda[k])
                          )
      i = i + 1
    }
  }
}

# vrednosti koje uzima slucajna velicina t_r
t_r_ocene = c()
for(i in 1:length(uzorci)) {
  R_ocena = mean(uzorci[[i]][1, ]) / mean(uzorci[[i]][2, ])
  t_r_ocene[i] = R_ocena * sum(br_razreda)
}
t_r_ocene

## [1] 888.8889 685.7143 900.0000 800.0000 977.7778 800.0000

Et_r = sum(t_r_ocene * 1/6)
Et_r

## [1] 842.0635

sum(br_petica) # nije nepristrasna - da bismo uporedjivali racunamo MSE

## [1] 850

MSE_t_r = sum((t_r_ocene - sum(br_petica))^2 / 6)
MSE_t_r

## [1] 8721.55

# ocena na osnovu psu bp
t_ocena_psu_bp = N * mean(uzorak_br_petica)
t_ocena_psu_bp

## [1] 600

# disperzija = MSE za psu bp
n = 2
D_t = N * (N-n) * var(br_petica) / n
D_t

## [1] 29166.67

MSE_t_r

## [1] 8721.55

# bolja kolicnicka ocena, ima manju srednje kvadratnu gresku

8. Дати су подаци за два стратума из популације величине 14:

I | x: 2, 2, 3, 4, 6, 7, 7, 3 | y: 6, 6, 7, 7, 8, 8, 6, 5

II | x: 10, 11, 12, 14, 15, 16 | y: 10, 11, 16, 17, 18, 18

Из првог стратума одабране су у узорак 4, а из другог 3 јединке. Вредности јединки одабраних у узорак су:

I | x: 2, 3, 7, 3 | y: 6, 7, 8, 5

II | x: 10, 14, 15 | y: 10, 17, 18

Наћи комбиновану и посебну количничку оцену средње вредности обележја популације, а затим одредити дисперзије тих оцена.

N = 14
stratumi = list()
stratumi[[1]] = rbind(c(2, 2, 3, 4, 6, 7, 7, 3), 
                      c(6, 6, 7, 7, 8, 8, 6, 5))
stratumi[[2]] = rbind(c(10, 11, 12, 14, 15, 16),
                      c(10, 11, 16, 17, 18, 18))

uzorak = list()
uzorak[[1]] = rbind(c(2, 3, 7, 3),
                    c(6, 7, 8, 5))
uzorak[[2]] = rbind(c(10, 14, 15),
                    c(10, 17, 18))

Ni = sapply(stratumi, ncol)
ni = sapply(uzorak, ncol)
xni = sapply(uzorak, function(x) {mean(x[1, ])})
yni = sapply(uzorak, function(x) {mean(x[2, ])})

xstr = sum(Ni*xni)
ystr = sum(Ni*yni)

Rc = xstr / ystr
Rc

## [1] 0.7605634

y = c(stratumi[[1]][2, ], stratumi[[2]][2, ])
xrc = Rc * mean(y)
xrc

## [1] 7.768612

x = c(stratumi[[1]][1, ], stratumi[[2]][1, ])

Ri = xni / yni
ty = sapply(stratumi, function(x) {sum(x[2, ])})
xrs = sum(Ri * ty) / N
xrs

## [1] 7.755495

si2x = sapply(stratumi, function(x) {var(x[1, ])})
si2y = sapply(stratumi, function(x) {var(x[2, ])})
rho = sapply(stratumi, function(x) {cor(x[1, ], x[2, ])})

R = mean(x)/mean(y)
Dxrc = sum(Ni^2*(1-ni/Ni) * (si2x + R^2*si2y - 2*R*rho*sqrt(si2x)*sqrt(si2y))/ni) / (N^2)
Dxrs = sum(Ni^2*(1-ni/Ni) * (si2x + Ri^2*si2y - 2*Ri*rho*sqrt(si2x)*sqrt(si2y))/ni) / (N^2)

Регресионо оцењивање

2. У квалификацијама за европско првенство у кошарци учествовало је 40 екипа. Одабран је прост случајан узорак од 6 екипа. Дати су подаци о броју постигнутих поена на свих 10 утакмица, као и број кошаркаша у тиму виших од 2 метра.

Број поена: 851, 945, 844, 978, 999, 821.

Број кошаркаша виших од 2m: 6, 4, 4, 6, 5, 7.

Познато је да је укупно на турниру било 211 играча виших од 2m. Регресионом методом оценити укупан број постигнутих поена у квалификацијама, ако је помоћно обележје број кошаркаша виших од 2 метра и ако је регресиони коефицијент:

а) Непознат, и оцењује се из узорка;

б) Познат, b0 = 20. У овом случају одредити оцену дисперзије оцене укупног броја поена.

N = 40
n = 6
br_poena = c(851, 945, 844, 978, 999, 821)
br_visih_od_2m = c(6, 4, 4, 6, 5, 7)

# suma pomocnog obelezja 
t_visih_od_2m = 211

rho_ocena = cor(br_poena, br_visih_od_2m)
rho_ocena

## [1] -0.2811843

b_ocena = rho_ocena * sqrt(var(br_poena))/sqrt(var(br_visih_od_2m))
b_ocena

## [1] -17.81818

# ili 
sum((br_poena - mean(br_poena)) * (br_visih_od_2m - mean(br_visih_od_2m))) / sum((br_visih_od_2m-mean(br_visih_od_2m))^2)

## [1] -17.81818

x_lr = mean(br_poena) + b_ocena * (t_visih_od_2m/N - mean(br_visih_od_2m))
t_lr = N * x_lr
t_lr

## [1] 36294.91

# b)
b0 = 20
x_lr_b0 = mean(br_poena) + b0 * (t_visih_od_2m/N - mean(br_visih_od_2m))
t_lr_b0 = N * x_lr_b0
t_lr_b0

## [1] 36206.67

D_x_lr_ocena = (1 - n/N) * (var(br_poena) + b0^2 * var(br_visih_od_2m) - 2 * b0 * rho_ocena * sqrt(var(br_poena)) * sqrt(var(br_visih_od_2m))) / n
D_x_lr_ocena

## [1] 1065.541

D_t_lr_ocena = N^2 * D_x_lr_ocena
D_t_lr_ocena

## [1] 1704866

5. Дата је популација од 12 јединки која је подељена на три стратума. Вредности обележја од интереса и помоћног обележја су:

I | x: 4, 3, 2, 1, 6 | y: 5, 5, 1, 2, 7

II | x: 2, 3 | y: 4, 5

III | x: 5, 3, 1, 7, 5 | y: 6, 3, 3, 2, 8

У узорак су одабрани елементи:

I | x: 4, 1, 6 | y: 5, 2, 7

II | x: 2 | y: 4

III | x: 5, 3 | y: 6, 3

а) Ако је познато да је b0 = 2, наћи комбиновану оцену средње вредности обележја x, а затим одредити дисперзију те оцене.

б) Ако је познато да је b10 = 2, b20 = 1, b30 = 3, наћи посебну оцену средње вредности обележја x, а затим одредити дисперзију те оцене и упоредити је са оценом из дела а).

N = 12

stratumi = list()
stratumi[[1]] = rbind(c(4, 3, 2, 1, 6),
                      c(5, 5, 1, 2, 7))
stratumi[[2]] = rbind(c(2, 3),
                      c(4, 5))
stratumi[[3]] = rbind(c(5, 3, 1, 7, 5),
                      c(6, 3, 3, 2, 8))
stratumi

## [[1]]
##      [,1] [,2] [,3] [,4] [,5]
## [1,]    4    3    2    1    6
## [2,]    5    5    1    2    7
## 
## [[2]]
##      [,1] [,2]
## [1,]    2    3
## [2,]    4    5
## 
## [[3]]
##      [,1] [,2] [,3] [,4] [,5]
## [1,]    5    3    1    7    5
## [2,]    6    3    3    2    8

uzorak = list()
uzorak[[1]] = rbind(c(4, 1, 6),
                    c(5, 2, 7))
uzorak[[2]] = rbind(c(2),
                    c(4))
uzorak[[3]] = rbind(c(5, 3),
                    c(6, 3))
uzorak

## [[1]]
##      [,1] [,2] [,3]
## [1,]    4    1    6
## [2,]    5    2    7
## 
## [[2]]
##      [,1]
## [1,]    2
## [2,]    4
## 
## [[3]]
##      [,1] [,2]
## [1,]    5    3
## [2,]    6    3

b0 = 2
Ni = sapply(stratumi, ncol)
ni = sapply(uzorak, ncol)

y_sr = mean(c(stratumi[[1]][2, ], stratumi[[2]][2, ],stratumi[[3]][2, ]))
y_sr

## [1] 4.25

xni = sapply(uzorak, function(x) {mean(x[1, ])})
yni = sapply(uzorak, function(x) {mean(x[2, ])})

x_str = sum(Ni * xni) / N
y_str = sum(Ni * yni) / N

x_lr_c = x_str + b0 * (y_sr - y_str)
x_lr_c

## [1] 3.055556

# vektor koji sadrzi koeficijente korelacije po stratumima
rho = sapply(stratumi, function(x) {cor(x[1, ],x[2, ])})
rho

## [1] 0.9020123 1.0000000 0.1572427

# vektor disperzija obelezja x po stratumima
s2_x = sapply(stratumi, function(x) {var(x[1, ])})
# vektor disperzija obelezja y po stratumima
s2_y = sapply(stratumi, function(x) {var(x[2, ])})

D_x_lr_c = sum(Ni^2 * (1 - ni/Ni) * (s2_x + b0^2 * s2_y - 2 * b0 * rho * sqrt(s2_x) * sqrt(s2_y)) / ni)
D_x_lr_c = D_x_lr_c / N^2
D_x_lr_c

## [1] 1.650463

# b)
b0_i = c(2, 1, 3)
y_sr_i = sapply(stratumi, function(x) {mean(x[2, ])})

x_lr_s = sum(Ni * (xni + b0_i * (y_sr_i - yni))) / N
x_lr_s

## [1] 2.930556

# srednja vrednost obelezja x na populaciji (vrednost koju ocenjujemo)
x_sr = mean(c(stratumi[[1]][1, ], stratumi[[2]][1, ], stratumi[[3]][1, ]))
x_sr

## [1] 3.5

x_lr_c

## [1] 3.055556

x_lr_s

## [1] 2.930556

D_x_lr_s = sum(Ni^2 * (1 - ni/Ni) * (s2_x + b0_i^2 * s2_y - 2 * b0_i * rho * sqrt(s2_x) * sqrt(s2_y)) / ni)
D_x_lr_s = D_x_lr_s / N^2
D_x_lr_s

## [1] 3.190394

D_x_lr_c < D_x_lr_s

## [1] TRUE

Кластер узорковање

1. Популација од 12 елемената подељена је на три кластера. Вредности обележја су:

I | 2, 4, 6

II | 3, 3, 4, 5, 5

III | 1, 1, 2, 4

Методом простог случајног узорковања без понављања у узорак су одабрани први и трећи кластер.Наћи непристрасну оцену суме обележја популације, дисперзију те оцене,а затим и количничку оцену суме обележја популације и дисперзију те оцене.

N = 3 # broj primarnih jedinica

klasteri = list()
klasteri[[1]] = c(2, 4, 6)
klasteri[[2]] = c(3, 3, 4, 5, 5)
klasteri[[3]] = c(1, 1, 2, 4)

# klasteri izabrani u uzorak
indeks = c(1, 3)

# sume obelezja po klasterima
ti = sapply(klasteri, sum)

n = 2
# ocena za t
t_ocena = N * sum(ti[indeks]) / n
t_ocena

## [1] 30

# vrednost koju smo ocenjivali
t = sum(ti)
t

## [1] 40

# disperzija ocene t_ocena
D_t_ocena = N^2 * (1 - n/N) * sum((ti - t/N)^2) / (n * (N - 1))
D_t_ocena

## [1] 56

# velicine klastera
Mi = sapply(klasteri, length)
Mu = sum(Mi)

# kolicnicka ocena
t_r = sum(ti[indeks]) * Mu / sum(Mi[indeks])
t_r

## [1] 34.28571

R = sum(ti) / sum(Mi)
D_t_r = N^2 * (1 - n/N) * sum((ti - R*Mi)^2) / (n * (N - 1))
D_t_r

## [1] 32.66667

4. Популација која садржи 250 елемената, издељена је на 20 скупина. Изабран је узорак од 5 скупина са вероватноћама пропорционалним величинама скупина (са понављањем), чије су величине и суме обележја дате:

Mi | 5 | 17 | 10 | 12 | 22

ti | 4 | 11 | 7 | 11 | 23

Наћи Хорвиц-Томпсонову оцену средње вредности обележја и оценити дисперзију те оцене.

Mu = 250
N = 20
n = 5

Mi = c(5, 17, 10, 12, 22)
ti = c(4, 11, 7, 11, 23)

# verovatnoce izbora u jednom biranju za klastere odabrane u uzorak
pi = Mi / Mu 

# verovatnoce ukljucenja 1.reda za klastere odabrane u uzorak
pii = 1 - (1 - pi)^n
pii

## [1] 0.0960792 0.2967989 0.1846273 0.2180396 0.3690801

# sve jedinke su razlicite
x_ht = sum(ti/pii) / Mu
x_ht

## [1] 0.9175012

D_x_ht_ocena = sum((1 - pii) * ti^2 / (pii^2))
for(i in 1:n) {
  for(j in 1:n) {
    if(i != j) {
      pij = pii[i] + pii[j] - 1 + (1 - pi[i] - pi[j])^n
      D_x_ht_ocena = D_x_ht_ocena + (pij - pii[i]*pii[j]) * (ti[i]*ti[j]) / (pii[i]*pii[j]*pij)
    }
  }
}
D_x_ht_ocena = D_x_ht_ocena / (Mu^2)
D_x_ht_ocena

## [1] 0.005581875

5. Дати су подаци о оценама из математике за 12 ученика који су подељени у 3 групе:

I | 2, 4, 5

II | 3, 3, 4, 5, 2

III| 1, 3, 5, 1

Бирањем пропорционалним величини групе, са понављањем, у узорак су одабране прва и трећа група. Наћи Хансен-Хурвицову оцену просечне оцене из математике и дисперзију те оцене.

klasteri = c()
klasteri[[1]] = c(2, 4, 5)
klasteri[[2]] = c(3, 3, 4, 5, 2)
klasteri[[3]] = c(1, 3, 5, 1)

Mi = c() # broj sekundarnih jedinica u klasteru
ti = c() # suma obelezja u klasteru
for(i in 1:length(klasteri)) {
  Mi[i] = length(klasteri[[i]])
  ti[i] = sum(klasteri[[i]])
}

Mu = 12
N = 3
n = 2
indeks = c(1, 3)

# HH ocena srednje vrednosti obelezja
x_hh = sum(ti[indeks] / Mi[indeks]) / n
x_hh

## [1] 3.083333

# prava vrednost
x_sr = mean(c(klasteri[[1]], klasteri[[2]], klasteri[[3]]))
x_sr

## [1] 3.166667

D_x_hh = sum(Mi * (ti/Mi - x_sr)^2) / (n*Mu)
D_x_hh

## [1] 0.1166667

Систематско узорковање

1. Популација се састоји од 20 јединки, чије су вредности обележја 1, 12, 4, 23, 23, 13, 15, 6, 16, 18, 35, 31, 23, 26, 19, 21, 32, 3, 27, 22. Формирати систематске узорке обима 5, а затим одредити дисперзију средине систематског узорка.

obelezje = c(1, 12, 4, 23, 23, 13, 15, 6, 16, 18, 35, 31, 23, 26, 19, 21, 32, 3, 27, 22)
N = length(obelezje)
n = 5

k = N / n
k

## [1] 4

# moguci uzorci obima n su dati po kolonama (ima ih k)
indeksi = matrix(1:N, n, k, byrow = T) 
indeksi

##      [,1] [,2] [,3] [,4]
## [1,]    1    2    3    4
## [2,]    5    6    7    8
## [3,]    9   10   11   12
## [4,]   13   14   15   16
## [5,]   17   18   19   20

uzorci = matrix(obelezje, n, k, byrow = T)
uzorci

##      [,1] [,2] [,3] [,4]
## [1,]    1   12    4   23
## [2,]   23   13   15    6
## [3,]   16   18   35   31
## [4,]   23   26   19   21
## [5,]   32    3   27   22

# disperzija
s2 = 0 
s_sis =  0
for(i in 1:k) {
  for (j in 1:n) {
    s2 = s2 + (uzorci[j, i] - mean(uzorci))^2
    s_sis = s_sis + (uzorci[j, i] - mean(uzorci[, i]))^2
  }
}
s2 = s2 / (N-1)
s2

## [1] 95.94737

s2 == var(obelezje)

## [1] TRUE

s_sis = s_sis / (k * (n - 1))
s_sis

## [1] 106.525

Dxsis = (N - 1) * s2 / N - k * (n - 1) * s_sis / N
Dxsis

## [1] 5.93

# ili 
xsr = mean(obelezje)
xisr = colMeans(uzorci)
mean((xisr - xsr)^2)

## [1] 5.93

2. Популација се састоји од 900 јединки. Одабрати узорак обима 200 користећи систематско узорковање.

N = 900
n = 200
N / n # nije ceo broj

## [1] 4.5

k = ceiling(N / n)
k

## [1] 5

# na slucajan nacin odaberemo jedan element od prvih k
r = sample(1:k, 1)
r # biramo u uzorak r, r + k, r + 2k, ...

## [1] 1

uzorak = seq(r, N, k) # uzimamo od r sa korakom k (najvise do N)
uzorak

##   [1]   1   6  11  16  21  26  31  36  41  46  51  56  61  66  71  76  81  86
##  [19]  91  96 101 106 111 116 121 126 131 136 141 146 151 156 161 166 171 176
##  [37] 181 186 191 196 201 206 211 216 221 226 231 236 241 246 251 256 261 266
##  [55] 271 276 281 286 291 296 301 306 311 316 321 326 331 336 341 346 351 356
##  [73] 361 366 371 376 381 386 391 396 401 406 411 416 421 426 431 436 441 446
##  [91] 451 456 461 466 471 476 481 486 491 496 501 506 511 516 521 526 531 536
## [109] 541 546 551 556 561 566 571 576 581 586 591 596 601 606 611 616 621 626
## [127] 631 636 641 646 651 656 661 666 671 676 681 686 691 696 701 706 711 716
## [145] 721 726 731 736 741 746 751 756 761 766 771 776 781 786 791 796 801 806
## [163] 811 816 821 826 831 836 841 846 851 856 861 866 871 876 881 886 891 896

length(uzorak)

## [1] 180

n - length(uzorak) # toliko nam fali

## [1] 20

max(uzorak)

## [1] 896

max(uzorak) + k # preko 900, ne moze u uzorak!

## [1] 901

# dopunjavamo do uzorka obima 200
dodatak = cumsum(c(max(uzorak) + k - N, rep(k, n - length(uzorak) - 1)))
dodatak

##  [1]  1  6 11 16 21 26 31 36 41 46 51 56 61 66 71 76 81 86 91 96

length(dodatak) == n - length(uzorak)

## [1] TRUE

uzorak = c(uzorak, dodatak)
uzorak # dopunjeni uzorak

##   [1]   1   6  11  16  21  26  31  36  41  46  51  56  61  66  71  76  81  86
##  [19]  91  96 101 106 111 116 121 126 131 136 141 146 151 156 161 166 171 176
##  [37] 181 186 191 196 201 206 211 216 221 226 231 236 241 246 251 256 261 266
##  [55] 271 276 281 286 291 296 301 306 311 316 321 326 331 336 341 346 351 356
##  [73] 361 366 371 376 381 386 391 396 401 406 411 416 421 426 431 436 441 446
##  [91] 451 456 461 466 471 476 481 486 491 496 501 506 511 516 521 526 531 536
## [109] 541 546 551 556 561 566 571 576 581 586 591 596 601 606 611 616 621 626
## [127] 631 636 641 646 651 656 661 666 671 676 681 686 691 696 701 706 711 716
## [145] 721 726 731 736 741 746 751 756 761 766 771 776 781 786 791 796 801 806
## [163] 811 816 821 826 831 836 841 846 851 856 861 866 871 876 881 886 891 896
## [181]   1   6  11  16  21  26  31  36  41  46  51  56  61  66  71  76  81  86
## [199]  91  96

length(uzorak)

## [1] 200