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

Називи расподела су: unif, norm, pois,beta, gamma, binom, geom, cauchy, chisq, t, exp, f, …

Одређени префикси се додају на име расподеле у зависности од тога шта желимо да израчунамо:

На пример за нормалну расподелу имамо: pnorm, dnorm, qnorm, rnorm.

Детаљније погледати на

# help("distribution")

Расподеле

Биномна расподела \(\mathcal{B}(n,p)\)

За рачунање појединачних вероватноћа: dbinom(x, size, prob, log = FALSE).

dbinom(0, size = 10, prob = 0.5)
## [1] 0.0009765625
dbinom(6, size = 10, prob = 0.5)
## [1] 0.2050781

Функција расподеле: pbinom(q, size, prob, lower.tail = T, log.p = F)

  • lower.tail подраѕумева да тражимо \(P\{X\le q\}\), да смо ставили FALSE, рачунала би се вероватноћа \(P\{X>q\}\).

  • log.p=T значило би да враћа \(log(p)\) уместо \(p\)

pbinom(8, size = 10, prob = 0.3)
## [1] 0.9998563

Функција квантила: qbinom(p, size, prob, lower.tail = T, log.p = F), где је \(q\), такво да је \(P\{X\le q\}=p\)

qbinom(0.99, size = 10, prob = 0.5)
## [1] 9

Узорак из биномне расподеле: rbinom(n, size, prob), где је \(n\) обим узорка.

rbinom(20, 10, 0.5)
##  [1] 5 3 1 5 4 5 7 2 6 5 5 6 8 3 3 5 5 5 3 7

График

x <- 0:50
plot(x, dbinom(x, size = 50, prob = 0.4))

plot(x, dbinom(x, size = 50, prob = 0.4), type = "h") #h od histogram

plot(x, dbinom(x, size = 50, prob = 0.4), type = "h", main = "Binomna raspodela B(50, 0.4)", xlab = "k", ylab = "binomne verovatnoce")

Биномни коефицијенти

choose(4, 2)
## [1] 6

Норамлна расподела \(\mathcal{N}(m,\sigma^2)\)

Густина: dnorm(x, mean = 0, sd = 1, log = F)

dnorm(3, mean = 2, sd = 4)
## [1] 0.09666703

Функција расподеле: pnorm(z, mean = 0, sd = 1, lower.tail = T, log.p = F)

pnorm(0)
## [1] 0.5
pnorm(3, mean = 4, lower.tail = F)
## [1] 0.8413447
pnorm(3, mean = 4, lower.tail = F, log.p = T)
## [1] -0.1727538

График густине

x <- seq(-5, 5, 0.01)
plot(x, dnorm(x), type = "l")

# Drugi nacin
curve(dnorm(x), from = -5, to = 5)

График функције расподеле

plot(x, pnorm(x), type = "l")

Квантили: qnorm(p, mean = 0, sd = 1)

qnorm(0.9)
## [1] 1.281552
qnorm(pnorm(2)) # Vraca bas 2, ovo su inverzne f-je 
## [1] 2

Узорак случајних бројева из \(\mathcal{N}(m,\sigma^2)\) расподеле: rnorm(size, mean = 0, sd = 1).

rnorm(20) #iz N(0,1)
##  [1] -1.49216064  0.12832243 -0.80986614  2.14963546 -1.12306738 -0.30941445
##  [7] -1.23402669  0.06910300 -0.41854419 -0.78307357  1.35472563  0.28097751
## [13]  0.14608470  0.02110867  1.27098269  2.16146584  1.00308353 -1.66557604
## [19] -0.62940704 -1.71222872

Униформна расподела \(\mathcal{U}[a,b]\)

Густина: dunif(x, min = а, max = b, log.p = F)

dunif(2)
## [1] 0
dunif(2, min = 3, max = 5)
## [1] 0

Функција расподеле: punif(q, min = 0, max = 1, lower.tail = TRUE, log.p = FALSE)

punif(1, min = 0, max = 3)
## [1] 0.3333333

Квантили: qunif(p, min = 0, max = 1,lower.tail = TRUE, log.p = FALSE)

qunif(0.5)
## [1] 0.5
runif(10) #iz U[0,1]
##  [1] 0.6568957 0.2093494 0.7565299 0.5215931 0.4359091 0.8630029 0.3976061
##  [8] 0.1734290 0.4003204 0.1214447

Експоненцијална расподела \(\mathcal{E}(\lambda)\)

Густина: dexp(x, rate = lambda, log = F)

Функција расподеле: pexp(q, rate = lambda, lower.tail = T, log.p = F)

Квантили: qexp(p, rate = lambda, lower.tail = T, log.p = F)

Случајни елемент: rexp(n, rate = lambda)

dexp(2, rate = 2)
## [1] 0.03663128
pexp(qexp(0.4))
## [1] 0.4
x <- seq(0, 25, 0.01)
y <- dexp(x, rate = 3)
plot(x, y, main = "Gustina E(3) raspodele", ylab = "f(x)", type = "l")

Студентова расподела

dt(x, df, ncp = 0, log = FALSE)

pt(q, df, ncp = 0, lower.tail = TRUE, log.p = FALSE)

qt(p, df, ncp = 0, lower.tail = TRUE, log.p = FALSE)

rt(n, df, ncp = 0)

dt(3, df = 3)
## [1] 0.02297204
rt(25, df = 6)
##  [1] -0.75324773 -1.07027740  0.09145127  0.45689804  0.48193766  1.48453876
##  [7]  3.77088105 -0.58145171 -0.54268695  0.78664026 -0.07828412  0.02593934
## [13]  0.21810974 -0.60667200 -1.12146541  3.71454012  1.55219229 -2.74984595
## [19] -0.79563309  0.51900210  0.14229703 -0.27378268  0.01161090 -0.58078160
## [25]  1.14755407

Хи-квадрат расподела

dchisq(x, df, ncp = 0, log = FALSE)

pchisq(q, df, ncp = 0, lower.tail = TRUE, log.p = FALSE)

qchisq(p, df, ncp = 0, lower.tail = TRUE, log.p = FALSE)

rchisq(n, df, ncp = 0)

Фишерова расподела

df(x, df1, df2, ncp, log = FALSE)

pf(q, df1, df2, ncp, lower.tail = TRUE, log.p = FALSE)

qf(p, df1, df2, ncp, lower.tail = TRUE, log.p = FALSE)

rf(n, df1, df2, ncp)

x <- seq(0, 100, 0.01)
curve(df(x, 12, 10), from = 0, to = 100) # crta gustinu

Пакет PASWR - Probability and Statistics with R

# install.packages("PASWR")
library(PASWR)
## Warning: package 'PASWR' was built under R version 4.1.3
## Loading required package: lattice

Неке функције у овом пакету везане су за вероватноћу.

На пример, функција bino.gen(samples, n, pi) за симулацију биномне расподеле.

bino.gen(10, 20, 0.4) # crta i teoretsku preko uzoracke raspodele

## $simulated.distribution
## Successes
##   9  10  11  12  13 
## 0.2 0.1 0.4 0.2 0.1 
## 
## $theoretical.distribution
##     0     1     2     3     4     5     6     7     8     9    10    11    12 
## 0.000 0.000 0.000 0.000 0.000 0.001 0.005 0.015 0.035 0.071 0.117 0.160 0.180 
##    13    14    15    16    17    18    19    20 
## 0.166 0.124 0.075 0.035 0.012 0.003 0.000 0.000

Функција Combinations(n, k) исписује све могуће комбинације.

Combinations(6, 2)
##   [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13] [,14]
##      1    1    2    1    2    3    1    2    3     4     1     2     3     4
## N    2    3    3    4    4    4    5    5    5     5     6     6     6     6
##   [,15]
##       5
## N     6

Функција SRS(Values, n) враћа све могуће узорке из популације Values обима n.

SRS(1:10, 5) #uzorak sa razlicitim elementima
##        [,1] [,2] [,3] [,4] [,5]
##   [1,]    1    2    3    4    5
##   [2,]    1    2    3    4    6
##   [3,]    1    2    3    4    7
##   [4,]    1    2    3    4    8
##   [5,]    1    2    3    4    9
##   [6,]    1    2    3    4   10
##   [7,]    1    2    3    5    6
##   [8,]    1    2    3    5    7
##   [9,]    1    2    3    5    8
##  [10,]    1    2    3    5    9
##  [11,]    1    2    3    5   10
##  [12,]    1    2    3    6    7
##  [13,]    1    2    3    6    8
##  [14,]    1    2    3    6    9
##  [15,]    1    2    3    6   10
##  [16,]    1    2    3    7    8
##  [17,]    1    2    3    7    9
##  [18,]    1    2    3    7   10
##  [19,]    1    2    3    8    9
##  [20,]    1    2    3    8   10
##  [21,]    1    2    3    9   10
##  [22,]    1    2    4    5    6
##  [23,]    1    2    4    5    7
##  [24,]    1    2    4    5    8
##  [25,]    1    2    4    5    9
##  [26,]    1    2    4    5   10
##  [27,]    1    2    4    6    7
##  [28,]    1    2    4    6    8
##  [29,]    1    2    4    6    9
##  [30,]    1    2    4    6   10
##  [31,]    1    2    4    7    8
##  [32,]    1    2    4    7    9
##  [33,]    1    2    4    7   10
##  [34,]    1    2    4    8    9
##  [35,]    1    2    4    8   10
##  [36,]    1    2    4    9   10
##  [37,]    1    2    5    6    7
##  [38,]    1    2    5    6    8
##  [39,]    1    2    5    6    9
##  [40,]    1    2    5    6   10
##  [41,]    1    2    5    7    8
##  [42,]    1    2    5    7    9
##  [43,]    1    2    5    7   10
##  [44,]    1    2    5    8    9
##  [45,]    1    2    5    8   10
##  [46,]    1    2    5    9   10
##  [47,]    1    2    6    7    8
##  [48,]    1    2    6    7    9
##  [49,]    1    2    6    7   10
##  [50,]    1    2    6    8    9
##  [51,]    1    2    6    8   10
##  [52,]    1    2    6    9   10
##  [53,]    1    2    7    8    9
##  [54,]    1    2    7    8   10
##  [55,]    1    2    7    9   10
##  [56,]    1    2    8    9   10
##  [57,]    1    3    4    5    6
##  [58,]    1    3    4    5    7
##  [59,]    1    3    4    5    8
##  [60,]    1    3    4    5    9
##  [61,]    1    3    4    5   10
##  [62,]    1    3    4    6    7
##  [63,]    1    3    4    6    8
##  [64,]    1    3    4    6    9
##  [65,]    1    3    4    6   10
##  [66,]    1    3    4    7    8
##  [67,]    1    3    4    7    9
##  [68,]    1    3    4    7   10
##  [69,]    1    3    4    8    9
##  [70,]    1    3    4    8   10
##  [71,]    1    3    4    9   10
##  [72,]    1    3    5    6    7
##  [73,]    1    3    5    6    8
##  [74,]    1    3    5    6    9
##  [75,]    1    3    5    6   10
##  [76,]    1    3    5    7    8
##  [77,]    1    3    5    7    9
##  [78,]    1    3    5    7   10
##  [79,]    1    3    5    8    9
##  [80,]    1    3    5    8   10
##  [81,]    1    3    5    9   10
##  [82,]    1    3    6    7    8
##  [83,]    1    3    6    7    9
##  [84,]    1    3    6    7   10
##  [85,]    1    3    6    8    9
##  [86,]    1    3    6    8   10
##  [87,]    1    3    6    9   10
##  [88,]    1    3    7    8    9
##  [89,]    1    3    7    8   10
##  [90,]    1    3    7    9   10
##  [91,]    1    3    8    9   10
##  [92,]    1    4    5    6    7
##  [93,]    1    4    5    6    8
##  [94,]    1    4    5    6    9
##  [95,]    1    4    5    6   10
##  [96,]    1    4    5    7    8
##  [97,]    1    4    5    7    9
##  [98,]    1    4    5    7   10
##  [99,]    1    4    5    8    9
## [100,]    1    4    5    8   10
## [101,]    1    4    5    9   10
## [102,]    1    4    6    7    8
## [103,]    1    4    6    7    9
## [104,]    1    4    6    7   10
## [105,]    1    4    6    8    9
## [106,]    1    4    6    8   10
## [107,]    1    4    6    9   10
## [108,]    1    4    7    8    9
## [109,]    1    4    7    8   10
## [110,]    1    4    7    9   10
## [111,]    1    4    8    9   10
## [112,]    1    5    6    7    8
## [113,]    1    5    6    7    9
## [114,]    1    5    6    7   10
## [115,]    1    5    6    8    9
## [116,]    1    5    6    8   10
## [117,]    1    5    6    9   10
## [118,]    1    5    7    8    9
## [119,]    1    5    7    8   10
## [120,]    1    5    7    9   10
## [121,]    1    5    8    9   10
## [122,]    1    6    7    8    9
## [123,]    1    6    7    8   10
## [124,]    1    6    7    9   10
## [125,]    1    6    8    9   10
## [126,]    1    7    8    9   10
## [127,]    2    3    4    5    6
## [128,]    2    3    4    5    7
## [129,]    2    3    4    5    8
## [130,]    2    3    4    5    9
## [131,]    2    3    4    5   10
## [132,]    2    3    4    6    7
## [133,]    2    3    4    6    8
## [134,]    2    3    4    6    9
## [135,]    2    3    4    6   10
## [136,]    2    3    4    7    8
## [137,]    2    3    4    7    9
## [138,]    2    3    4    7   10
## [139,]    2    3    4    8    9
## [140,]    2    3    4    8   10
## [141,]    2    3    4    9   10
## [142,]    2    3    5    6    7
## [143,]    2    3    5    6    8
## [144,]    2    3    5    6    9
## [145,]    2    3    5    6   10
## [146,]    2    3    5    7    8
## [147,]    2    3    5    7    9
## [148,]    2    3    5    7   10
## [149,]    2    3    5    8    9
## [150,]    2    3    5    8   10
## [151,]    2    3    5    9   10
## [152,]    2    3    6    7    8
## [153,]    2    3    6    7    9
## [154,]    2    3    6    7   10
## [155,]    2    3    6    8    9
## [156,]    2    3    6    8   10
## [157,]    2    3    6    9   10
## [158,]    2    3    7    8    9
## [159,]    2    3    7    8   10
## [160,]    2    3    7    9   10
## [161,]    2    3    8    9   10
## [162,]    2    4    5    6    7
## [163,]    2    4    5    6    8
## [164,]    2    4    5    6    9
## [165,]    2    4    5    6   10
## [166,]    2    4    5    7    8
## [167,]    2    4    5    7    9
## [168,]    2    4    5    7   10
## [169,]    2    4    5    8    9
## [170,]    2    4    5    8   10
## [171,]    2    4    5    9   10
## [172,]    2    4    6    7    8
## [173,]    2    4    6    7    9
## [174,]    2    4    6    7   10
## [175,]    2    4    6    8    9
## [176,]    2    4    6    8   10
## [177,]    2    4    6    9   10
## [178,]    2    4    7    8    9
## [179,]    2    4    7    8   10
## [180,]    2    4    7    9   10
## [181,]    2    4    8    9   10
## [182,]    2    5    6    7    8
## [183,]    2    5    6    7    9
## [184,]    2    5    6    7   10
## [185,]    2    5    6    8    9
## [186,]    2    5    6    8   10
## [187,]    2    5    6    9   10
## [188,]    2    5    7    8    9
## [189,]    2    5    7    8   10
## [190,]    2    5    7    9   10
## [191,]    2    5    8    9   10
## [192,]    2    6    7    8    9
## [193,]    2    6    7    8   10
## [194,]    2    6    7    9   10
## [195,]    2    6    8    9   10
## [196,]    2    7    8    9   10
## [197,]    3    4    5    6    7
## [198,]    3    4    5    6    8
## [199,]    3    4    5    6    9
## [200,]    3    4    5    6   10
## [201,]    3    4    5    7    8
## [202,]    3    4    5    7    9
## [203,]    3    4    5    7   10
## [204,]    3    4    5    8    9
## [205,]    3    4    5    8   10
## [206,]    3    4    5    9   10
## [207,]    3    4    6    7    8
## [208,]    3    4    6    7    9
## [209,]    3    4    6    7   10
## [210,]    3    4    6    8    9
## [211,]    3    4    6    8   10
## [212,]    3    4    6    9   10
## [213,]    3    4    7    8    9
## [214,]    3    4    7    8   10
## [215,]    3    4    7    9   10
## [216,]    3    4    8    9   10
## [217,]    3    5    6    7    8
## [218,]    3    5    6    7    9
## [219,]    3    5    6    7   10
## [220,]    3    5    6    8    9
## [221,]    3    5    6    8   10
## [222,]    3    5    6    9   10
## [223,]    3    5    7    8    9
## [224,]    3    5    7    8   10
## [225,]    3    5    7    9   10
## [226,]    3    5    8    9   10
## [227,]    3    6    7    8    9
## [228,]    3    6    7    8   10
## [229,]    3    6    7    9   10
## [230,]    3    6    8    9   10
## [231,]    3    7    8    9   10
## [232,]    4    5    6    7    8
## [233,]    4    5    6    7    9
## [234,]    4    5    6    7   10
## [235,]    4    5    6    8    9
## [236,]    4    5    6    8   10
## [237,]    4    5    6    9   10
## [238,]    4    5    7    8    9
## [239,]    4    5    7    8   10
## [240,]    4    5    7    9   10
## [241,]    4    5    8    9   10
## [242,]    4    6    7    8    9
## [243,]    4    6    7    8   10
## [244,]    4    6    7    9   10
## [245,]    4    6    8    9   10
## [246,]    4    7    8    9   10
## [247,]    5    6    7    8    9
## [248,]    5    6    7    8   10
## [249,]    5    6    7    9   10
## [250,]    5    6    8    9   10
## [251,]    5    7    8    9   10
## [252,]    6    7    8    9   10

Основне статистичке функције

x <- rnorm(50)
mean(x) # Uzoracka sredina
## [1] -0.06783661
median(x) # Uzoracka medijana
## [1] -0.181151
sd(x) # Standardna devijacija
## [1] 0.9696448
var(x) # Popravljena uzoracka disperzija
## [1] 0.9402111
quantile(x) # Funkcija koja vraca kvantile (0%, 25%, 50%, 75%, 100%)
##         0%        25%        50%        75%       100% 
## -1.9949972 -0.7877706 -0.1811510  0.6950214  2.3525893
p <- seq(0, 1, 0.1) # Ako zelimo druge verovatnoce, prosledjujemo ih u vektoru
quantile(x, p)
##         0%        10%        20%        30%        40%        50%        60% 
## -1.9949972 -1.1619705 -0.9780135 -0.7009013 -0.3723898 -0.1811510  0.1268643 
##        70%        80%        90%       100% 
##  0.3902563  0.7626456  1.0251396  2.3525893
range(x) # Raspon uzorka
## [1] -1.994997  2.352589

Подаци са недостајућим вредностима

attach(airquality)
mean(Ozone)# ne moze da izracuna zbog NA vrednosti, pa zato
## [1] NA
mean(Ozone, na.rm = TRUE)
## [1] 42.12931
summary(Ozone) # vidimo da ima 321 NA vrednost
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
##    1.00   18.00   31.50   42.13   63.25  168.00      37

Графички приказ података

Хистограм

Хистограм је графичка репрезентација низа датих нумеричких података. Користи се за оцену густине.

Хистограм се формира на следећи начин:

  1. Добијени подаци се сортирају.

  2. Одабере се дужина подеока \(d\).

  3. Подели се цео интервал (распон података) на подинтервале дужине \(d\).

  4. На \(x\)-оси се означе ти добијени интервали, а одговарајућа вредност на \(y\)-оси је број елемената из узорка који се налазе у том интервалу.

У \(R\)-у постоји уграђена функција hist(). Аргументи ове функције су:

  • x - вектор чији хистограм желимо да прикажемо

  • breaks - вектор са крајевима сваког интервала.

  • freq - да ли су апсолутне или релативне фреквенције, у другом случају укупно површина је 1

  • right - ако је TRUE ћелије су облика \((a,b]\), а за FALSE обрнуто

Такође, постоје и аргументи col, border, main, xlab, ylab, xlim, ylim.

Кад правимо хистограм прво треба да одаберемо колико ћемо категорија (подинтервала) да имамо. Тај број се добија из формуле \[k=[log_2(N)]+1,\]

где је \(k\) број категорија, а \(N\) обим узорка, тј. величина тог вектора.

Одавде можемо наћи ширину сваког интервала по формули \(d=\frac{R}{k}\), где је \(R\) распон узорка (разлика између највећег и најмањег члана).

  • Пример: Тврди се да је просечна цена безоловног бензина у Америци била \(1.35\$\). У рекламне сврхе компанија жели да покаже како је њихова цена нижа. Да би поткрепили своју тврдњу, статистичари из фирме су сакупили следеће податке на основу случајног узорка:
(cene <- c(1.22, 1.37, 1.27, 1.20, 1.42, 1.41, 1.22, 1.24, 1.28, 1.42, 1.48, 1.32, 1.40, 1.26, 1.39, 1.45, 1.44, 1.49, 1.47, 1.47, 1.24, 1.34, 1.27, 1.35, 1.34, 1.45, 1.49, 1.45, 1.23, 1.20, 1.42, 1.34,   1.43, 1.21, 1.49, 1.36, 1.24, 1.20, 1.45,  1.23, 1.25, 1.24, 1.35, 1.23, 1.39, 1.38, 1.46, 1.48, 1.26, 1.36, 1.22, 1.46, 1.39,  1.22, 1.29, 1.47, 1.24, 1.35, 1.21, 1.21))
##  [1] 1.22 1.37 1.27 1.20 1.42 1.41 1.22 1.24 1.28 1.42 1.48 1.32 1.40 1.26 1.39
## [16] 1.45 1.44 1.49 1.47 1.47 1.24 1.34 1.27 1.35 1.34 1.45 1.49 1.45 1.23 1.20
## [31] 1.42 1.34 1.43 1.21 1.49 1.36 1.24 1.20 1.45 1.23 1.25 1.24 1.35 1.23 1.39
## [46] 1.38 1.46 1.48 1.26 1.36 1.22 1.46 1.39 1.22 1.29 1.47 1.24 1.35 1.21 1.21

Израчунати у \(R\)-у узорачку средину и медијану, узорачки распон и нацртати хистограм над прикупљеним подацима.

mean(cene) 
## [1] 1.340167
median(cene) 
## [1] 1.35
range(cene) # ili c(min(cene),max(cene))
## [1] 1.20 1.49

Правимо хистограм

n <- length(cene)
k <- floor(log(n, base = 2)) + 1
d <- diff(range(cene)) / k
k
## [1] 6
d
## [1] 0.04833333
cene <- sort(cene) # Sortiramo vektor

Правимо поделу на интервале

podela <- cene[1] + 0:k * d
podela
## [1] 1.200000 1.248333 1.296667 1.345000 1.393333 1.441667 1.490000
hist(cene,breaks = podela, col = 'coral', main = "Histogram")

Упоредимо овај хистограм са оним који се добија ако не задамо сами поделе.

par(mfrow = c(1, 2))
hist(cene, breaks = podela, col = 'coral', main = "Histogram cena")
hist(cene, col = 'lightblue' , main = "Histogram cena")

Полигон фреквенција

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

Да бисмо конструисали полигон фреквенција, треба на почетку да поделимо податке на неке класе (интервале), баш као код хистограма.

Потом се означи тачка на графику чија \(x\) координата има вредност средине интервала, а \(y\) координата је фреквенција одговарајуће класе. Полигон добијамо тако што дужима спајамо те означене тачке. Можемо укључити једну класу пре најмање вредности међу подацима и једну после највеће. На тај начин добијамо да полигон додирује \(x\) осу са обе стране.

Показаћемо полигон на примеру података из базе discoveries која садржи податке о броју великих изума и научних открића између 1860. и 1959. године.

discoveries
## Time Series:
## Start = 1860 
## End = 1959 
## Frequency = 1 
##   [1]  5  3  0  2  0  3  2  3  6  1  2  1  2  1  3  3  3  5  2  4  4  0  2  3  7
##  [26] 12  3 10  9  2  3  7  7  2  3  3  6  2  4  3  5  2  2  4  0  4  2  5  2  3
##  [51]  3  6  5  8  3  6  6  0  5  2  2  2  6  3  4  4  2  2  4  7  5  3  3  0  2
##  [76]  2  2  1  3  4  2  2  1  1  1  2  1  4  4  3  2  1  4  1  1  1  0  0  2  0

По формули добијамо да треба да имамо отприлике 7 категорија, па за ширину интервала можемо узети 2 и додајемо још две класе на крајевима.

range(discoveries)
## [1]  0 12
cut.points <- seq(-2, 14, by = 2)

Табеле фреквенција

disc.cut <- cut(discoveries, cut.points, right = FALSE)
disc.freq <- table(disc.cut)
cbind(disc.freq)
##         disc.freq
## [-2,0)          0
## [0,2)          21
## [2,4)          46
## [4,6)          19
## [6,8)          10
## [8,10)          2
## [10,12)         1
## [12,14)         1

Цртамо хистограм и полигон фреквенција

hist(discoveries, 
     breaks = cut.points, 
     col = "slategray3", 
     border = "dodgerblue4",
     right = FALSE,
     xlab = "x-osa", 
     main = "Histogram")

plot(disc.freq,
     type = "b",
     col = "orange",
     main = "Poligon frekvencija")

Хистограм као оцена густине

Неколико примера у којима ћемо да генеришемо узорак из познате расподеле и да упоредимо облик хистограма са теоријском густином. (Користимо хистограм из пакета lattice).

library (lattice)
  • Узорак обима 500 из стандардне нормалне расподеле
x <- rnorm(500)
hist (x, prob = TRUE, xlab = "", ylab = "", col ='coral1', border = 'bisque',  main = "N(0,1)")
curve(dnorm(x), add = TRUE, lwd = 3, col = 'cornflowerblue')

  • Узорак обима 1000 из експоненцијалне расподеле са параметром 4
y <- rexp(1000, 4)
hist (y, prob = TRUE, xlab = "", ylab = "", col = 'cornflowerblue', border = 'bisque',  main = "Exp(1)")
curve(dexp(x, 4), add = TRUE, lwd = 3, col = 'coral1')

  • Узорак из униформне \(\mathcal{U}(0 ,1)\) расподеле обима 500
z <- runif (500)
hist(z, prob = TRUE, xlab = "", ylab = "", col = 'darksalmon', border  = 'bisque', main = "U(0,1)")
curve(dunif(x), add = TRUE, lwd = 3, col = 'cadetblue')

  • Пример: Нека је \(S_n=X_1+...+X_n\), где су \(X_i,\) \(i=1,...,n\) независне случајне величине са експоненцијалном \(E(2)\) расподелом. Генерисати случајну величину \(S_n\). Генерисати 1000 бројева из исте расподеле као \(S_n\) и нацртати њихов хистограм. Да ли је резултат очекиван? (Подсетите се ЦГТ)
s_n <- function(n) {
  x <- rexp(n, rate = 2)
  sum(x)
}
N <- 1000
n <- 100
s <- replicate(N, s_n(n))
hist(s,probability = T, col = 'lightblue', main = "")

# Znamo da je EX=1/lambda, Dx=1/lambda^2
EX <- 1/2
DX <- 1/4

Хоћемо да стандардизујемо податке

s.z <- (s-n*EX)/sqrt(n*DX)
hist(s.z, probability = T, col = 'lightblue', main = "")
curve(dnorm(x), add = T, lwd = 2, col = 'coral')

Специјални случај централне граничне теореме је Муавр-Лапласова теорема:

Познато нам је да се биномна случајна величина \(\mathcal{B}(n,p)\) може представити као збир \(n\) независних случајних величина са Бернулијевом \(\mathcal{Ber}(p)\) расподелом (збир \(n\) независних индикатора) и такође \(ES_n=np\), \(DS_n=np(1−p)\). Ако је \(np>10\), можемо је апроксимирати нормалном.

Генерисаћемо \(\mathcal{B}(n,p)\) као збир \(n\) независних Бернулијевих случајних величина.

binom <- function(n, p) {
  x <- sample(c(0,1), n, replace = TRUE, prob = c(1-p, p))
  b <- sum(x)
  return(b)
}
uzorak <- replicate(1000, binom(50, 0.7))
# standardizujemo podatke
st.uzorak <- (uzorak - 50 * 0.7) / sqrt(50 * 0.7 * 0.3)
plot(density(st.uzorak), lwd = 2, col = 'coral', main = "")  
curve(dnorm(x), col = 'lightblue', lwd = 2, add = T, from = -4, to = 4)

Хистограм за груписане податке

library(ISwR)
## Warning: package 'ISwR' was built under R version 4.1.3
data(energy)
attach(energy)
head(energy)
##   expend stature
## 1   9.21   obese
## 2   7.53    lean
## 3   7.48    lean
## 4   8.08    lean
## 5   8.09    lean
## 6  10.15    lean
# Delimo podatke u dve grupe
expend.lean <- expend[stature == "lean"]
expend.obese <- expend[stature == "obese"]
hist(expend.lean) # ovde breaks prosledjujemo kao broj celija koje hocemo, ali R zadrzava pravo da ga izmeni da bi granice bile sto zaokruzeniji brojevi

hist(expend.obese)

Box-plot

Користи се за детекцију аутлајера.

boxplot(discoveries, col = "pink")

range(discoveries, na.rm = T)
## [1]  0 12
median(discoveries, na.rm = T)   
## [1] 3
q1 <- quantile(discoveries, na.rm = T)[2]
q3 <- quantile(discoveries, na.rm = T)[4]
qr <- IQR(discoveries, na.rm = T)
q1
## 25% 
##   2
q3
## 75% 
##   4
q1 - 1.5 * qr # f1
## 25% 
##  -1
q3 + 1.5 * qr # f3
## 75% 
##   7
q1 - 3 * qr # F1
## 25% 
##  -4
q3 + 3 * qr # F3
## 75% 
##  10
max(discoveries[discoveries < q3 + 1.5 * qr], na.rm = T) # a3
## [1] 6

Boxplot за груписане податке

boxplot(expend ~ stature) # y~x, za neke vektore x i y znaci: y je opisano preko x, i zove se model formula

# drugi nacin
boxplot(expend.lean, expend.obese) # samo ih shvata kao odvojene vektore, pa crta i odvojene plotove

Емпиријска функција расподеле

Нека је \((X_1,...,X_n)\) узорак из расподеле \(F\). Функција \[F_n(x)=\frac{1}{n}\sum_{k=1}^nI\{X_k\le x\}\] је емпиријска функција расподеле.

Детаљније погледати у пакету stepfun.

Пример

x <- rpois(100, lambda = 3)
table(x)
## x
##  0  1  2  3  4  5  6  7  8 11 
##  6 16 23 24 16  6  4  2  2  1
Fn <- ecdf(x)
plot(Fn, main = "Empirijska funkcija raspodele")
curve(ppois(x, lambda = 3), from = 0, to = 10, col = "red", add = TRUE)

x <- runif(100)
Fn <- ecdf(x)
plot(Fn, main = "Empirijska funkcija raspodele")
curve(punif(x), from = -1, to = 2, add = T, col = "blue")

Задатак

Користећи емпиријску функцију расподеле испитати да ли дати подаци имају експоненцијалну расподелу.

Y <- c(56, 83, 104, 116, 244, 305, 429, 452, 453, 503, 552, 614, 661, 673, 683, 685, 753, 763, 806, 834, 838, 862, 897, 904, 981, 1007, 1008, 1049, 1060, 1107, 1125, 1141, 1153, 1154, 1193, 1201, 1253, 1313, 1329, 1347, 1454, 1464, 1490, 1491, 1532, 1549, 1568, 1574, 1586, 1599, 1608, 1723, 1769, 1795, 1927, 1957, 2005, 2010, 2016, 2022, 2037, 2065, 2096, 2139, 2150, 2156, 2160, 2190, 2210, 2220, 2248, 2285, 2325, 2337, 2351, 2437, 2454, 2546, 2565, 2584, 2624, 2675, 2701, 2755, 2877, 2879, 2922, 2986, 3092, 3160, 3185, 3191, 3439, 3617, 3685,3756, 3826, 3995, 4007, 4159, 4300, 4487, 5074, 5579, 5623, 6869, 7739)
plot(ecdf(Y), main = "Empirijska funkcija raspodele")
curve(pexp(x, rate = 1/mean(Y)), from = 0, to = 8000, add = T, col = "blue")

Q-Q plot (quantile versus quantile)

График који се добија када се \(k\)-та вредност по величини плотује са очекиваном \(k\)-том вредношћу по величини у нормалној расподели. Ако се добије приближно права линија, јесте нормална расподела.

qqnorm(x)
qqline(x, col = "red")

set.seed(5432)
x1 <- rnorm(10000)
qqnorm(x1)
qqline(x1, col = "red")

x2 <- rnorm(10000, 5, 16)
qqnorm(x2)
qqline(x2, col = "red")