I’ll use an example where \({X1, X2, X3, X4}\) will be the number of independent random variables, \(n = 4\), uniformly distributed on integers from 1 to 6, \(k=6\).
To calculate the number of ways \(n\) is uniformly distributed on \(k\), we’ll do \(6^4\), which gives us 1,296. This number will be the denominator, \(k^n\).
When \(Y=1\), the number of combinations where \(X1\) to \(X4\) are greater than 1 is \(5^4\). The total number of combinations is \(6^4\), so therefore to find the distribution when \(Y=1\) is \(6^4 - 5^4\), which can be expressed as \(k^n - (k-1)^n\).
When \(Y=2\), the number of combinations where \(X1\) to \(X4\) are greater than 2 is \(4^4\). Using the methods from the previous calculations, we can express the distribution when \(Y=2\) as \(6^4 - (6^4 - 5^4) - 4^4\). This can be expressed as \(k^n - (k^n - ((k-1)^n) - (k-2)^n\), simplified as \((k-1)^n - (k-2)^n\) when \(Y=2\).
Where \(P(Y=y)\), we can therefore state that when \(y=2\), \((k-2+1)^n - (k-2)^n)\), with \((k-2+1)^n\) simplified as \((k-1)^n\). The final distribution of \(Y\) is:
\((k-y+1)^n - (k-y)^n) / k^n\)
Reference: https://math.dartmouth.edu/archive/m20f10/public_html/HW5Solutions.pdf
# Probability of Failure in 10 years
p <- 1/10
# Geometric Model for machine failing after 8 years
geom_model <- 1 - pgeom(8, p)
geom_model
## [1] 0.3874205
Expected Value:
geom_expected_value <- 1/p
geom_expected_value
## [1] 10
Standard Deviation:
sd_geom <- sqrt(1-p)/p
sd_geom
## [1] 9.486833
# rate of machine failing after 10 years
lambda <- 1/10
exp_model <- 1 - pexp(8, lambda)
exp_model
## [1] 0.449329
Expected Value:
exp_expected_value <- 1/lambda
exp_expected_value
## [1] 10
Standard Deviation:
sd_exp <- sqrt(1/lambda^2)
sd_exp
## [1] 10
# Number of trials (using number of years machine will fail)
n <- 8
# Probability of 1 failure after 8 years
binom_model <- choose(n, 0) * p^0 * (1-p)^n
binom_model
## [1] 0.4304672
Expected Value:
# p = 0.1
binom_expected_value <- n * p
binom_expected_value
## [1] 0.8
Standard Deviation:
sd_binom <- sqrt(binom_expected_value*(1-p))
sd_binom
## [1] 0.8485281
n <- 8
p <- 0.1
lambda_rate <- n * p
pois_model <- ppois(0,lambda_rate)
pois_model
## [1] 0.449329
Expected Value:
pois_expected_value <- lambda_rate
pois_expected_value
## [1] 0.8
Standard Deviation:
sd_pois <- sqrt(lambda_rate)
sd_pois
## [1] 0.8944272