# Chapter 4: Probability Distributions¶

## Example, page-83¶

In [1]:
#Variable Declaration
from scipy import *
from pylab import *

x1 = array([1,2,3,4])            # Domain of f(x)
x2 = array([0,1,2,3,4])          # Domain of h(x)
f = (x1-2)/2.0                   # list of values of f(x)
h = (x2*x2)/25.0                 # list of values of h(x)

#Calculation

#Results
if(sum(f)!=1):
print "f(x) can not serve as probability distribution"

if(sum(h)!=1):
print "h(x) can not serve as probability distribution"

print "f(0) is negative(P>=0) so f(x) can not serve as probability distribution"

h(x) can not serve as probability distribution
f(0) is negative(P>=0) so f(x) can not serve as probability distribution


## Example, page-87¶

In [2]:
#Variable Declaration
n = 5            # total number of installations
p = 0.6          # probability of heat installations reduced by one-third

#Calculation
from scipy import *
from pylab import *

def fact(n):
"""returns factorial of number n"""
if(n==1 or n==0):
return 1
else:
return n*fact(n-1)

def comb(n,r):
"""returns number of total combination of selecting "r" items out of "n" """
return fact(n)/(fact(r)*fact(n-r))
p1 = comb(5,4) * (p**4) * ((1-p)**1)                                                # probability of exactly four out of five
p1=round(p1,3)

p2 = (comb(5,4) * (p**4) * ((1-p)**1)) + (comb(5,5) * (p**5) * ((1-p)**0))          # probability of atleast four out of five
p2=round(p2,3)

#Results
print "probability of reduction of utility bill in exactly 4 out of 5: ",p1
print "probability of reduction of utility bill in atleast 4 out of 5: ",p2

probability of reduction of utility bill in exactly 4 out of 5:  0.259
probability of reduction of utility bill in atleast 4 out of 5:  0.337


## Example, page-87¶

In [3]:
#Variable Declaration
p = 0.05          # probability of failure due to load
n = 16            # total columns

#Calculation
from scipy import *
from pylab import *
def fact(n):
"""returns factorial of number n"""
if(n==1 or n==0):
return 1
else:
return n*fact(n-1)

def comb(n,r):
"""returns number of total combination of selecting "r" items out of "n" """
return fact(n)/(fact(r)*fact(n-r))
p1 = comb(16,0)*(p**0)*((1-p)**16) + comb(16,1)*(p**1)*((1-p)**15) + comb(16,2)*(p**2)*((1-p)**14)   # probability of failure of at most 2
p1 = round(p1,3)
p2 = comb(16,0)*(p**0)*((1-p)**16) + comb(16,1)*(p**1)*((1-p)**15) + comb(16,2)*(p**2)*((1-p)**14) + comb(16,3)*(p**3)*((1-p)**13)     # probability of at most 3
p2 = 1-p2                   # probability of failures of at least 4
p2 = round(p2,3)

#Results
print "probability of failure of at most 2 out of 16 columns: ",p1
print "probability of failure of at least 4 out of 16 columns: ",p2

probability of failure of at most 2 out of 16 columns:  0.957
probability of failure of at least 4 out of 16 columns:  0.007


## Example, Page-88¶

In [4]:
#Variable Declaration
p = 0.3         # probability of not passing inspection
n = 18          # total panels

#Calculation
from scipy import *
from pylab import *

def fact(n):
"""returns factorial of number n"""
if(n==1 or n==0):
return 1
else:
return n*fact(n-1)

def comb(n,r):
"""returns number of total combination of selecting "r" items out of "n" """
return fact(n)/(fact(r)*fact(n-r))

p1 = comb(18,6)*(p**6)*((1-p)**12)
p1 = round(p1,4)

#Results
print "Required probability: ",p1

Required probability:  0.1873


## Example, Page-88¶

In [5]:
#Variable Declaration
p = 0.1         # probability for repairing
n = 20          # total hard drives
p1 = 0

#Calculation
from scipy import *
from pylab import *

def fact(n):
"""returns factorial of number n"""
if(n==1 or n==0):
return 1
else:
return n*fact(n-1)

def comb(n,r):
"""returns number of total combination of selecting "r" items out of "n" """
return fact(n)/(fact(r)*fact(n-r))

for i in range(0,5):
p1 = p1 + comb(20,i)*(p**i)*((1-p)**(n-i))
p1 = 1 - round(p1,4)

#Results
print "As the probability(",p1,") is very less, so reject the claim."

As the probability( 0.0432 ) is very less, so reject the claim.


## Example, page-90¶

In [6]:
#Variable Declaration
n = 20    # total chargers
n1 = 10   # total selected chargers
r = 5     # total defective chargers
r1 = 2    # defective chargers to be selected out of 5

#Calculation
from scipy import *
from pylab import *

def fact(n):
"""returns factorial of number n"""
if(n==1 or n==0):
return 1
else:
return n*fact(n-1)

def comb(n,r):
"""returns number of total combination of selecting "r" items out of "n" """
return fact(n)/(fact(r)*fact(n-r))
p = (comb(r,r1)*comb(n-r,n1-r1)) / float(comb(n,n1))               # Required probability
p = round(p,3)

#Results
print "probability of exactly 2 defective out of 10: ",p

probability of exactly 2 defective out of 10:  0.348


## Example, Page-90¶

In [7]:
#Variable Declaration
n = 100    # total chargers
n1 = 10   # total selected chargers
r = 25     # total defective chargers
r1 = 2    # defective chargers to be selected out of 5
p = 0.25

#Calculation
from scipy import *
from pylab import *

def fact(n):
"""returns factorial of number n"""
if(n==1 or n==0):
return 1
else:
return n*fact(n-1)

def comb(n,r):
"""returns number of total combination of selecting "r" items out of "n" """
return fact(n)/(fact(r)*fact(n-r))
p1 = (comb(r,r1)*comb(n-r,n1-r1)) / float(comb(n,n1))
p1 = round(p1,3)

p2 = comb(n1,r1)*(p**r1)*((1-p)**(n1-r1))
p2 = round(p2,3)

#Results
print "For Hypergeometric distribution: ",p1
print "For Binomial distribution: ",p2

For Hypergeometric distribution:  0.292
For Binomial distribution:  0.282


## Example, page-95¶

In [8]:
#Variable Declaration
n = 3   # no of tosses
p0 = (0.5)*(0.5)*(0.5)       # probability of Number of heads=0
p1 = 3*(0.5)*(0.5)*(0.5)     # probability of Number of heads=1
p2 = 3*(0.5)*(0.5)*(0.5)     # probability of Number of heads=2
p3 = (0.5)*(0.5)*(0.5)       # probability of Number of heads=3

#Calculation
from scipy import *
from pylab import *
# As we know, mean = Sum(x*P(x))
mean_head = (p0)*0 + (p1)*1 + (p2)*2 + (p3)*3    # mean of number of heads

#Results

Mean of number of heads in 3 flips:  1.5


## Example, Page-95¶

In [9]:
#Variable Declaration
a = ([[0,0.18], [1,0.5], [2,0.29], [3,0.03]])
U = 0

#Calculation
from scipy import *
from pylab import *

for i in range(0,4):
U = U + a[i][0]*a[i][1]

#Results
print "Mean of the probability distribution: ",U

Mean of the probability distribution:  1.17


## Example, Page-96¶

In [10]:
#Variable Declaration
n = 3    # No. of flips
p = 0.5  # probability of head in one flip

#Calculation
from scipy import *
from pylab import *

#Results
print "Mean Number of heads: ",n*p

Mean Number of heads:  1.5


## Example, page-96¶

In [11]:
#Variable Declaration
N = 20       # Total car chargers
n = 10       # Number of selected car chargers
a = 5        # Defective car chargers

#Calculation
from scipy import *
from pylab import *
# As we know, for hypergeometric distribution, mean = (n*a)/N
Mean=float((n*a))/N           # Mean of probability distribution of number of defectives

#Results
print "mean of probability distribution of the number of defectives in a sample of 10 car chargers: ",Mean

mean of probability distribution of the number of defectives in a sample of 10 car chargers:  2.5


## Example, page-97¶

In [12]:
#Variable declaration
n1 = 4        # for subplot-1
p1 = 0.5      # for subplot-1
n2 = 16       # for subplot-2
p2 = 0.5      # for subplot-2

#Calculation
from scipy import *
from pylab import *
from scipy.integrate import *

%matplotlib inline
def fact(n):
"""returns factorial of number n"""
if(n==1 or n==0):
return 1
else:
return n*fact(n-1)

def comb(n,r):
"""returns number of total combination of selecting "r" items out of "n" """
return fact(n)/(fact(r)*fact(n-r))
# To calculate standard deviation of subplot-1

x1 = [0,1,2,3,4]                 # for subplot-1
f1=[]
Mean1 = (n1)*(p1)                # mean of subplot-1
y1=0

for each in x1:
y1 = y1 + ((each-Mean1)**2)*(comb(n1,each)*((p1)**n1))

std_dev1 = sqrt(y1)              # standard deviation of subplot-1

for each in x1:
f1.append((comb(n1,each)*((p1)**n1)))

# To calculate standard deviation of subplot-2

x2=[2,3,4,5,6,7,8,9,10,11,12,13,14]      # for subplot-2
f2=[]
Mean2 = (n2)*(p2)                        # Mean of subplot-2
y2 = 0

for each in x2:
y2 = y2 + ((each-Mean2)**2)*(comb(n2,each)*((p2)**n2))

std_dev2 = sqrt(y2)              # Standard deviation of subplot-2
std_dev2=round(std_dev2,1)

for each in x2:
f2.append((comb(n2,each)*((p2)**n2)))

# subplot-1
subplot(2,1,1)
bar(x1,f1,width=1.0,fill=False)
title("Probability distribution plot-1")
xlabel("$x$")
ylabel("$f1(x)$")

# subplot-2
subplot(2,1,2)
title("Probability distribution plot-2")
bar(x2,f2,width=1.0,fill=False)
xlabel("$x$")
ylabel("$f2(x)$")

#Results
print "First subplot standard deviation: ",std_dev1,"\nSecond subplot standard deviation: ",std_dev2

First subplot standard deviation:  1.0
Second subplot standard deviation:  2.0


## Example, Page-98¶

In [13]:
#Variable Declaration
n = 16
p = 0.5

#Calculation
from scipy import *
from pylab import *

Variance = float(n*p*(1-p))                # Variance of Sample
std_dev=sqrt(Variance)                     # standard deviation of sample

#Results
print "standard deviation of given sample: ",int(std_dev)

standard deviation of given sample:  2


## Example, page-98¶

In [14]:
#Variable Declaration
N = 20       # Total car chargers
n = 10       # Number of selected car chargers
a = 5        # Defective car chargers

#Calculation
from scipy import *
from pylab import *
# As we know, for hypergeometric distribution, Variance = (n.a/n)*(1-a/N)*((N-n)/(N-1))
Variance = (float((n*a))/N)*(1-(float(a)/N))*((float(N-n)/(N-1)))                # Variance of Sample
std_dev=sqrt(Variance)                                                           # standard deviation of sample
std_dev = round(std_dev,2)

#Results
print "standard deviation of given sample: ",std_dev

standard deviation of given sample:  0.99


## Example, page-99¶

In [15]:
#Variable Declaration
U = (1.0/6)*(1+2+3+4+5+6)                                  # Mean of probability distribution
U1 = (1.0/6)*(1*1 + 2*2 + 3*3 + 4*4 + 5*5 + 6*6)           # Second moment about the origin

#Calculation
from scipy import *
from pylab import *
# As we know, Variance = (U1)-U*U

Variance = (U1) - U*U                                      # Variance of probability distribution
Variance = round(Variance,2)

#Results
print "Variance of Given probability distribution: ",Variance,"or 35/12"

Variance of Given probability distribution:  2.92 or 35/12


## Example, page-99¶

In [16]:
#Variable declaration

#Calculation
from scipy import *
from pylab import *
from scipy.integrate import *

l = array([[0,0.05],[1,0.20],[2,0.45],[3,0.20],[4,0.10]])
x = l[:,0]
f = l[:,1]

U = sum(x*f)                # Mean of distribution
U2 = sum(x*x*f)             # Second moment about origin
Variance = U2 - (U)*(U)     # Variance of given probability distribution

Variance = round(Variance,2)
U = round(U,2)

#Results
print "Mean of given probability distribution: ",U
print "Variance of given probability distribution: ",Variance

Mean of given probability distribution:  2.1
Variance of given probability distribution:  0.99


## Example, page-101¶

In [17]:
#Variable Declaration
# as we know K = (x-Mean)/std_dev,       P(|x-Mean| < K*std_dev) >= 1-(1/(K*K))

K1 = (28-18)/2.5               # Corresponding to 28 customers
K2 = (18-8)/2.5                # Corresponding to 8 customers

#Calculation
from scipy import *
from pylab import *
# P( 8 < X < 28) >= 1-(1/(K*K))

K = 4      # K=K1=K2
P = 1-(1.0 / (K*K))            # Required probability

#Results
print " Probability of number of customers greater than 8 but less than 28 is greater than or equal to : ",P,"or 15/16"

 Probability of number of customers greater than 8 but less than 28 is greater than or equal to :  0.9375 or 15/16


## Example, Page-102¶

In [28]:
#Variable Declaration
n = 40000
p = 0.5

#Calculation
from scipy import *
from pylab import *

Mean = n*p
Variance = float(n*p*(1-p))                # Variance of Sample
std_dev=sqrt(Variance)                     # standard deviation of sample

k_square = 1/(1-0.99)
k = sqrt(k_square)

#Results
print "Probability limit: (",round(p1,3),",",round(p2,3),")"

Probability limit: ( 0.475 , 0.525 )


## Example, page-105¶

In [19]:
#Variable Declaration
r = 2           # defective books
n = 100         # total books
p = 0.05        # probability of a book being defective

#Calculation
from scipy import *
from pylab import *
def fact(n):
"""returns factorial of number n"""
if(n==1 or n==0):
return 1
else:
return n*fact(n-1)

def comb(n,r):
"""returns number of total combination of selecting "r" items out of "n" """
return fact(n)/(fact(r)*fact(n-r))
# using bionomial distribution
p1 = comb(100,2)*pow(p,2)*pow(1-p,98)                         # required probability using bionomial distribution
p1=round(p1,3)
# using poisson distribution
lamda = n*p
# we know P = (pow(e,-lamda)*pow(lamda,r)) / fact(r)
p2 = (pow(e,-lamda)*pow(lamda,r)) / fact(r)                 # required probability
p2 = round(p2,3)

#Results
print "using bionomial distribution: probability = ",p1
print "using poisson distribution: probability = ",p2

using bionomial distribution: probability =  0.081
using poisson distribution: probability =  0.084


## Example, Page-106¶

In [27]:
#Variable Declaration
n = 3840
p = 1.0/1200
l = []
x = [0,1,2,3,4,5,6,7,8,9,10]

#Calculation
from scipy import *
from pylab import *
%matplotlib inline

def fact(n):
"""returns factorial of number n"""
if(n==1 or n==0):
return 1
else:
return n*fact(n-1)

lamda = n*p

for i in range(0,11):
val = (pow(e,-lamda)*pow(lamda,i))/fact(i)
l.append(val)

#Results
print "Lamda: ",lamda
bar(x,l,width=1,fill=False)
xlabel("$No. of Failed generators$")
ylabel("$Probability$")

Lamda:  3.2

Out[27]:
<matplotlib.text.Text at 0xa350898>

## Example, Page-106¶

In [21]:
#Variable Declaration

#Calculation
from scipy import *
from pylab import *

def fact(n):
"""returns factorial of number n"""
if(n==1 or n==0):
return 1
else:
return n*fact(n-1)

lamda = 1.3

val = (pow(e,-lamda)*pow(lamda,0))/fact(0)
p = 1 - val
variance = lamda

#Results
print "Lamda: ",lamda
print "Probability(X>=1): ",round(p,3)
print "Variance: ",variance

Lamda:  1.3
Probability(X>=1):  0.727
Variance:  1.3


## Example, page-107¶

In [22]:
#Variable Declaration
alpha = 6    # no of bad checks per day

#Calculation
from scipy import *
from pylab import *
def fact(n):
"""returns factorial of number n"""
if(n==1 or n==0):
return 1
else:
return n*fact(n-1)
# 4 bad checks on any day
lamda = alpha*1                                        # As T=1

p1 = (pow(e,-lamda)*pow(lamda,4)) / fact(4)           # required probability
p1 = round(p1,3)

# 10 bad checks over any 2 consecutive days
lamda = alpha*2

p2 = ((pow(e,-lamda)*pow(lamda,10)) / fact(10))       # required probability
p2 = round(p2,3)

#Results
print " probability of 4 bad cheacks on any day: ",p1
print " probability of 10 bad cheacks on any 2 consecutive days: ",p2

 probability of 4 bad cheacks on any day:  0.134
probability of 10 bad cheacks on any 2 consecutive days:  0.105


## Example, Page-108¶

In [23]:
#Variable Declaration
l = []

#Calculation
from scipy import *
from pylab import *

def fact(n):
"""returns factorial of number n"""
if(n==1 or n==0):
return 1
else:
return n*fact(n-1)

lamda = (0.2)*3
l.append((pow(e,-lamda)*pow(lamda,1))/fact(1))

lamda = (0.2)*5
val = (pow(e,-lamda)*pow(lamda,0))/fact(0) + (pow(e,-lamda)*pow(lamda,1))/fact(1)
l.append(1 - val)

lamda = (0.2)*15
l.append((pow(e,-lamda)*pow(lamda,0))/fact(0) + (pow(e,-lamda)*pow(lamda,1))/fact(1))

#Results
print "Part(a): ",round(l[0],3)
print "Part(b): ",round(l[1],3)
print "Part(c): ",round(l[2],3)

Part(a):  0.329
Part(b):  0.264
Part(c):  0.199


## Example, page-109¶

In [24]:
#Variable Declaration
p = 0.05                         # probability of any device to show excessive drift

#Calculation
from scipy import *
from pylab import *
p1 = pow(1-p,5)*pow(p,1)         # required probability
p1 = round(p1,3)

#Results
print "probability of sixth measuring device to be the first to show excessive drift: ",p1

probability of sixth measuring device to be the first to show excessive drift:  0.039


## Example, Page-112¶

In [25]:
#Variable Declaration
n = 8
n1 = 2
n2 = 5
n3 = 1
p1 = 0.3
p2 = 0.5
p3 = 0.2

#Calculation
from scipy import *
from pylab import *

def fact(n):
"""returns factorial of number n"""
if(n==1 or n==0):
return 1
else:
return n*fact(n-1)

p = (fact(8)/(fact(2)*fact(5)*fact(1)))*(pow(p1,2)*pow(p2,5)*pow(p3,1))

#Results
print "Required Probability: ",round(p,4)

Required Probability:  0.0945