Factorials is writen as \(n!\) , and calculated by multiplied all integers from one to \(n\).
kaijou <- function(x) #first of all, define the function
{
if(x == round(x)) #check integer or not. "round(x)" shows the nearest integer to x.
{
if(x >= 1){
y <- c(1:x) # make vector from one to x. By doing this, we can use function "prod(x)".
z <- prod(y) #multiplies all numbers in the vector x by using "prod(x)"
return(z)
}
if(x == 0){ # if the number is zero, above code can't work well. So I need to make another code for the case the number is zero.
a <- 1
return(1)
}
if(x < 1){
print("n is inappropriate number!") # if the number is under zero, display about that
}
}
else{
print("n must be integer!") # if the number is not integer, display about that
}
}
The formula to calculate conbination is
So,write code along this formula.
kumiawase <- function(n,k) #define function. n must be larger than k
{
if(n>=k)
{
m <- kaijou(n)/{kaijou(n-k)*kaijou(k)}
return(m)
}
if(n < k){ #if n is smaller than k, show that
print("n must be larger than k!!")
}
}
To calculate the PMF of binominal distributions, we use the formula this.
\(f(x)=\begin{pmatrix} n \\ x \end{pmatrix}p^x(1-p)^{n-x}\)
Along this formula, calculate the PMF of binomial distribution by using fanction of conbinations that made in Q2.
nikoubunpu <- function(x,y,z){
if( 1>= z && z >= 0 && y >= x) #z is probability, so z must be between zero and one. Besides, y must not be under x.
{
t <- kumiawase(y,x)*z^x*(1-z)^(y-x)
return(t)
}
else{
print("the formula is unsuitable!!") #If y is under x, show that.
}
}
n <- 2
probability <- 0.5
f1 <- rep(NA, n+1) #The number of trials are "n + 1" times, because we calculate the PM from "x = 0" to "x = n"
b1 <- 0:n
for (i in 1:(n+1)) { #the range of i starts from one, not zero. So the end of the range of i must be "n + 1" to align the ranges.
p <- nikoubunpu(b1[i], n, probability) #b1[i] shows i-th number between zero and n, I defined b1 four lines on this sentence to do this.
f1[i] <- p
}
#This roop repeats the calcuculation between x = 0, the first number and x = n,(n + 1)th number.
plot(x = b1, y = f1, xlim = c(0,n), ylim = c(0,max(f1)), type = "h", lwd = 2, xlab = "x", ylab = "The Value of PM")
title(main = " the PM of binomial distributions for n=2, p=0.5")
n <- 2
probability <- 0.3
f1 <- rep(NA, n+1)
b1 <- 0:n
for (i in 1:(n+1)) {
p <- nikoubunpu(b1[i], n, probability)
f1[i] <- p
}
plot(x = b1, y = f1, xlim = c(0,n), ylim = c(0,max(f1)), type = "h", lwd = 2, xlab = "x", ylab = "The Value of PM")
title(main = " the PM of binomial distributions for n=2, p=0.3")
n <- 5
probability <- 0.5
f1 <- rep(NA, n+1)
b1 <- 0:n
for (i in 1:(n+1)) {
p <- nikoubunpu(b1[i], n, probability)
f1[i] <- p
}
plot(x = b1, y = f1, xlim = c(0,n), ylim = c(0,max(f1)), type = "h", lwd = 2, xlab = "x", ylab = "The Value of PM")
title(main = " the PM of binomial distributions for n=5, p=0.5")
n <- 5
probability <- 0.3
f1 <- rep(NA, n+1)
b1 <- 0:n
for (i in 1:(n+1)) {
p <- nikoubunpu(b1[i], n, probability)
f1[i] <- p
}
plot(x = b1, y = f1, xlim = c(0,n), ylim = c(0,max(f1)), type = "h", lwd = 2, xlab = "x", ylab = "The Value of PM")
title(main = " the PM of binomial distributions for n=5, p=0.3")
n <- 10
probability <- 0.5
f1 <- rep(NA, n+1)
b1 <- 0:n
for (i in 1:(n+1)) {
p <- nikoubunpu(b1[i], n, probability)
f1[i] <- p
}
plot(x = b1, y = f1, xlim = c(0,n), ylim = c(0,max(f1)), type = "h", lwd = 2, xlab = "x", ylab = "The Value of PM")
title(main = " the PM of binomial distributions for n=10, p=0.5")
n <- 10
probability <- 0.3
f1 <- rep(NA, n+1)
b1 <- 0:n
for (i in 1:(n+1)) {
p <- nikoubunpu(b1[i], n, probability)
f1[i] <- p
}
plot(x = b1, y = f1, xlim = c(0,n), ylim = c(0,max(f1)), type = "h", lwd = 2, xlab = "x", ylab = "The Value of PM")
title(main = " the PM of binomial distributions for n=10, p=0.3")
n <- 50
probability <- 0.5
f1 <- rep(NA, n+1)
b1 <- 0:n
for (i in 1:(n+1)) {
p <- nikoubunpu(b1[i], n, probability)
f1[i] <- p
}
plot(x = b1, y = f1, xlim = c(0,n), ylim = c(0,max(f1)), type = "h", lwd = 2, xlab = "x", ylab = "The Value of PM")
title(main = " the PM of binomial distributions for n=50, p=0.5")
n <- 50
probability <- 0.3
f1 <- rep(NA, n+1)
b1 <- 0:n
for (i in 1:(n+1)) {
p <- nikoubunpu(b1[i], n, probability)
f1[i] <- p
}
plot(x = b1, y = f1, xlim = c(0,n), ylim = c(0,max(f1)), type = "h", lwd = 2, xlab = "x", ylab = "The Value of PM")
title(main = " the PM of binomial distributions for n=50, p=0.3")