2116.sagews

import scipy.stats

Generating 10 values of the normal distribution with mean 2 and standard deviation 3.

scipy.stats.norm.rvs(2,3,size=10)
sum(scipy.stats.norm.rvs(2,3,size=10))

Traceback (most recent call last): File "<stdin>", line 1, in <module> File "/home/sage/sagenb/sage_notebook/worksheets/jason3/215/code/4.py", line 8, in <module> exec compile(ur'mean(scipy.stats.norm.rvs(_sage_const_2 ,_sage_const_3 ,size=_sage_const_10 ))' + '\n', '', 'single') File "", line 1, in <module> NameError: name 'mean' is not defined

r.rnorm(10,2,3).sage()

[-0.85415578501983902, 2.27374894289966, 0.52189865593636198, -1.49242045832038, 5.9421766167807997, -1.33639158610723, -4.1697116302391199, 1.4192608556292801, 2.4232202358326602, -1.98620114112697]

dist=RealDistribution('gaussian',3)
[dist.get_random_element()+2 for _ in range(10)]

[4.46748691304, 4.04462061546, 1.00466636762, -0.971824030927, -4.03142330473, -1.19382494784, 8.28026107492, 4.08986379227, 2.47245674343, -1.66389184829]

finance.TimeSeries(10).randomize('normal',2,3)

[5.8754, 3.7232, 5.4661, 5.8238, -5.2359, 5.9938, 1.1857, 6.7851, 4.9460, 2.2808]

Generating 10 values of the uniform distribution from 100 to 1000.  For scipy, we specify the starting value and the range.  For R, we specify the starting and ending values.

scipy.stats.uniform.rvs(100,900,size=10)

array([ 773.55044344, 666.58526535, 215.09909016, 527.02901491, 289.8517689 , 337.44003005, 790.19643141, 864.45231476, 196.01508513, 991.30528085])

r.runif(10,100,1000).sage()

[979.26270824391395, 916.69925677124399, 126.272868458182, 940.17426606733397, 102.23680264316501, 418.27806057408498, 552.32447043526895, 935.21712596993905, 419.11423190031201, 458.37259285617603]

dist=RealDistribution('uniform',[100,1000])
[dist.get_random_element() for _ in range(10)]

[970.326854521, 910.559309227, 592.50902792, 255.42579249, 975.415918441, 870.058850385, 743.334395043, 648.132038303, 727.955941902, 637.800585129]

finance.TimeSeries(10).randomize('uniform',100,900)

[866.1378, 277.1640, 361.6241, 238.2037, 132.3171, 380.5158, 887.4178, 987.7804, 768.9125, 594.8248]

Generating 10 values of the Bernoulli distribution with one of the values having probability 0.5 (i.e., a coin-flip).  A success is a 1, a failure is a 0.

scipy.stats.bernoulli.rvs(0.5,size=10)

array([0, 0, 0, 1, 0, 1, 1, 1, 1, 0])

Generating 10 values of the Binomial distribution.  The values below model 30 coin tosses, with the probability of success of each trial being 0.5.  The values represent the number of successful (head) flips in the 30 tosses.

scipy.stats.binom.rvs(30,0.5,size=10)

array([16, 15, 11, 13, 13, 18, 13, 14, 11, 17])

r.rbinom(10,30,0.5).sage()

[18, 14, 15, 15, 12, 11, 15, 15, 17, 18]

Comparing distributions

bin=scipy.stats.binom(30, 0.5)

bin.rvs(100)

array([15, 14, 16, 12, 15, 18, 15, 13, 17, 16, 18, 16, 19, 12, 17, 13, 18, 16, 13, 18, 15, 15, 15, 20, 14, 16, 15, 20, 14, 15, 13, 15, 11, 11, 15, 21, 18, 23, 19, 14, 15, 16, 14, 16, 20, 14, 19, 19, 11, 19, 20, 17, 15, 14, 18, 12, 16, 10, 13, 14, 16, 16, 13, 16, 17, 15, 18, 15, 5, 21, 16, 12, 11, 18, 18, 18, 17, 14, 11, 19, 14, 16, 17, 15, 13, 15, 13, 16, 16, 19, 15, 15, 18, 14, 13, 16, 16, 12, 15, 17])

normal_dist=scipy.stats.norm.rvs(50,25,size=1000)
uniform_dist=scipy.stats.uniform.rvs(0,100,size=1000)
binom_dist=scipy.stats.binom.rvs(100,0.5,size=1000)

import matplotlib.pyplot as plt

plt.figure()
plt.boxplot([normal_dist, uniform_dist, binom_dist])
plt.title("Boxplot")
plt.xlabel("Distribution")
plt.ylabel("Values")
plt.savefig('test.png')

plot(scipy.stats.norm(0,1).pdf,(-10,10))+plot(scipy.stats.uniform(-5,10).pdf,(-10,10))

scipy.stats.norm(0,1)

<scipy.stats.distributions.rv_frozen object at 0x4e9d390>

Suppose we have a 5% failure rate in producing widgets.  What is a histogram of the probability of having $k$ failures in a box of 20 widgets?

binom_dist=scipy.stats.binom(20,.05) # 20 trials, "success" (the failure of the bulb) 5% of the time.

binom_dist.pmf(2) # the probability of 2 successes (i.e., 2 bulb failures)

0.18867680126765396

bar_chart([binom_dist.pmf(x) for x in range(21)])

a=[binom_dist.pmf(x) for x in range(21)]

a

[0.35848592240854188, 0.37735360253530764, 0.18867680126765396, 0.05958214776873294, 0.013327585685111276, 0.002244646010123974, 0.00029534815922693802, 3.1089279918572466e-05, 2.6589515720321089e-06, 1.8659309275470548e-07, 1.0802758021455361e-08, 5.168783179243519e-10, 2.0403123635048814e-11, 6.6080474425689317e-13, 1.7430501486614958e-14, 2.2204460492503131e-16, 0.0, 0.0, 0.0, 0.0, 1.1102230246251565e-16]

sum(a)/len(a)

0.047619047619047616