Simulating Bernoulli(p) Distribution
p = 0.2;
n = 50;
U = runif(n, min = 0, max = 1);
X = sum(U < p);
X
## [1] 10
p1=0.35
p2=0.15
p3=0.4
p4=0.1
P = c(p1,p1+p2, p1+p2+p3, p1+p2+p3+p4)
X = c(1, 2, 3, 4)
counter = 1
set.seed(2021)
r = runif(1, min = 0, max = 1);
while(r > P[counter])
counter = counter + 1
end
## function (x, ...)
## UseMethod("end")
## <bytecode: 0x00000253752021b8>
## <environment: namespace:stats>
X[counter]
## [1] 2
U = runif(1, min = 0, max = 1);
lambda = 1
X = -(1/lambda)*log(U)
Interpretation:
Normal distribution
dnorm(0)
## [1] 0.3989423
dnorm(1)
## [1] 0.2419707
dnorm(0, mean=4, sd=3)
## [1] 0.05467002
pnorm(0)
## [1] 0.5
pnorm(1)
## [1] 0.8413447
x = seq(-3,3, by=0.1); y = pnorm(x) ;plot(x,y, type="l")
pnorm(0, lower.tail=FALSE)
## [1] 0.5
pnorm(1, lower.tail=FALSE)
## [1] 0.1586553
x = seq(-3,3, by=0.1); y = pnorm(x, lower.tail=FALSE)
plot(x,y, type="l")
qnorm(0.68); qnorm(0.95);qnorm(0.997)
## [1] 0.4676988
## [1] 1.644854
## [1] 2.747781
x = seq(0,1, by=0.05); y = qnorm(x);plot(x,y, type="l")
x=rnorm(1000)
hist(x)
pnorm(1) - pnorm(-1) # within one standard deviation
## [1] 0.6826895
pnorm(2) - pnorm(-2) # within two standard deviation
## [1] 0.9544997
pnorm(3) - pnorm(-3) # within three standard deviation
## [1] 0.9973002
Are the sum of rolls normal ?
DiceR = read.csv("Dice.csv", header=T)
names(DiceR)= "Sum"
DiceR$Sum = as.numeric(DiceR$Sum)
## Warning: NAs introduced by coercion
Dice = na.omit(DiceR$Sum)
hist(Dice)
summary(Dice)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 7.00 15.00 17.50 17.23 20.00 27.00
sd(Dice)
## [1] 3.670899
cs = (Dice-mean(Dice))/(sd(Dice))
mean(cs)
## [1] 2.261974e-16
sd(cs)
## [1] 1
onesdcs = cs[cs >-1& cs <1] # within one sd
twosdcs = cs[cs >-2& cs <2] # within two sd
threesdcs = cs[cs >-3& cs <3]# within three sd
length(onesdcs)/length(cs)
## [1] 0.6653846
length(twosdcs)/length(cs)
## [1] 0.9576923
length(threesdcs)/length(cs)
## [1] 1
library(moments)
Kurtosis = kurtosis(Dice)
Kurtosis
## [1] 3.002322
Skewness = skewness(Dice)
Skewness
## [1] -0.2548971
normaldata = replicate(1000, rnorm(97,0,1),
simplify=FALSE)
library(moments)
mean(sapply(normaldata, kurtosis))
## [1] 2.947374
mean(sapply(normaldata, skewness))
## [1] 0.0006454242
normaldata = replicate(1000, rnorm(97,0,1),
simplify=FALSE)
library(moments)
hist(sapply(normaldata, kurtosis))
x=rnorm(1000)
qqnorm(x)
qqline(x)
qqnorm(Dice)
qqline(Dice)