Exercise from Think Stats, 2nd Edition (thinkstats2.com)
Allen Downey

Read the pregnancy file.


In [3]:
%matplotlib inline
import matplotlib
matplotlib.style.use('ggplot')
import nsfg
preg = nsfg.ReadFemPreg()
preg


nsfg.py:42: SettingWithCopyWarning: 
A value is trying to be set on a copy of a slice from a DataFrame

See the the caveats in the documentation: http://pandas.pydata.org/pandas-docs/stable/indexing.html#indexing-view-versus-copy
  df.birthwgt_lb[df.birthwgt_lb > 20] = np.nan
Out[3]:
caseid pregordr howpreg_n howpreg_p moscurrp nowprgdk pregend1 pregend2 nbrnaliv multbrth ... laborfor_i religion_i metro_i basewgt adj_mod_basewgt finalwgt secu_p sest cmintvw totalwgt_lb
0 1 1 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 3410.389399 3869.349602 6448.271112 2 9 NaN 8.8125
1 1 2 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 3410.389399 3869.349602 6448.271112 2 9 NaN 7.8750
2 2 1 NaN NaN NaN NaN 5 NaN 3 5 ... 0 0 0 7226.301740 8567.549110 12999.542264 2 12 NaN 9.1250
3 2 2 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 7226.301740 8567.549110 12999.542264 2 12 NaN 7.0000
4 2 3 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 7226.301740 8567.549110 12999.542264 2 12 NaN 6.1875
5 6 1 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 4870.926435 5325.196999 8874.440799 1 23 NaN 8.5625
6 6 2 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 4870.926435 5325.196999 8874.440799 1 23 NaN 9.5625
7 6 3 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 4870.926435 5325.196999 8874.440799 1 23 NaN 8.3750
8 7 1 NaN NaN NaN NaN 5 NaN 1 NaN ... 0 0 0 3409.579565 3787.539000 6911.879921 2 14 NaN 7.5625
9 7 2 NaN NaN NaN NaN 5 NaN 1 NaN ... 0 0 0 3409.579565 3787.539000 6911.879921 2 14 NaN 6.6250
10 12 1 NaN NaN NaN NaN 5 NaN 1 NaN ... 0 0 0 3612.781968 4146.013572 6909.331618 1 31 NaN 7.8125
11 14 1 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 2418.069494 2810.302771 3039.904507 2 56 NaN 7.0000
12 14 2 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 2418.069494 2810.302771 3039.904507 2 56 NaN 4.0000
13 14 3 NaN NaN NaN NaN 3 NaN NaN NaN ... 0 0 0 2418.069494 2810.302771 3039.904507 2 56 NaN NaN
14 15 1 NaN NaN NaN NaN 1 NaN NaN NaN ... 0 0 0 1667.816099 3200.862017 5553.495599 1 33 NaN NaN
15 15 2 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 1667.816099 3200.862017 5553.495599 1 33 NaN 7.6875
16 15 3 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 1667.816099 3200.862017 5553.495599 1 33 NaN 7.5000
17 18 1 NaN NaN NaN NaN 5 NaN 1 NaN ... 0 0 0 2957.257457 3404.403067 4153.371741 2 14 NaN 6.3125
18 18 2 NaN NaN NaN NaN 1 NaN NaN NaN ... 0 0 0 2957.257457 3404.403067 4153.371741 2 14 NaN NaN
19 21 1 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 3408.342437 3965.763949 7237.122630 1 48 NaN 8.7500
20 21 2 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 3408.342437 3965.763949 7237.122630 1 48 NaN 8.1875
21 23 1 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 6210.373020 8120.841310 13533.382043 2 64 NaN 5.5625
22 23 2 NaN NaN NaN NaN 1 NaN NaN NaN ... 0 0 0 6210.373020 8120.841310 13533.382043 2 64 NaN NaN
23 24 1 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 3409.573258 4068.628645 7424.840414 1 27 NaN 6.7500
24 24 2 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 3409.573258 4068.628645 7424.840414 1 27 NaN 7.3750
25 24 3 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 3409.573258 4068.628645 7424.840414 1 27 NaN 6.8125
26 28 1 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 3407.794208 3808.343516 6949.846082 2 57 NaN 8.1250
27 31 1 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 3405.679025 4272.084519 5211.943113 1 2 NaN 7.1250
28 31 2 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 3405.679025 4272.084519 5211.943113 1 2 NaN 6.0625
29 31 3 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 3405.679025 4272.084519 5211.943113 1 2 NaN 7.4375
... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ...
13563 12547 2 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 3453.545517 6628.022524 11499.619080 1 52 NaN 7.6875
13564 12547 3 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 3453.545517 6628.022524 11499.619080 1 52 NaN 7.6250
13565 12550 1 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 3080.452699 3745.326058 5268.550165 1 79 NaN 8.1250
13566 12551 1 NaN NaN NaN NaN 5 NaN 1 NaN ... 0 0 0 2418.538866 3653.453268 3951.940400 2 75 NaN 7.5000
13567 12554 1 NaN NaN NaN NaN 3 NaN NaN NaN ... 0 0 0 1914.676604 2177.957240 2764.045534 2 75 NaN NaN
13568 12554 2 NaN NaN NaN NaN 4 NaN NaN NaN ... 0 0 0 1914.676604 2177.957240 2764.045534 2 75 NaN NaN
13569 12556 1 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 2474.619764 3250.573384 3965.699528 1 44 NaN 5.8125
13570 12556 2 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 2474.619764 3250.573384 3965.699528 1 44 NaN 6.6875
13571 12556 3 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 2474.619764 3250.573384 3965.699528 1 44 NaN 6.0000
13572 12556 4 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 2474.619764 3250.573384 3965.699528 1 44 NaN 5.8125
13573 12561 1 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 2418.089703 2698.650781 4497.301527 1 10 NaN 6.5625
13574 12561 2 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 2418.089703 2698.650781 4497.301527 1 10 NaN 6.1250
13575 12564 1 NaN NaN NaN NaN 3 NaN NaN NaN ... 0 0 0 1820.850938 2129.214067 2768.191208 2 44 NaN NaN
13576 12565 1 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 3195.641221 3834.241709 6652.409365 1 78 NaN 6.4375
13577 12565 2 35 1 8 NaN NaN NaN NaN NaN ... 0 0 0 3195.641221 3834.241709 6652.409365 1 78 NaN NaN
13578 12566 1 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 2080.317155 2422.820274 2627.548587 2 2 NaN 6.0000
13579 12566 2 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 2080.317155 2422.820274 2627.548587 2 2 NaN 7.0000
13580 12568 1 NaN NaN NaN NaN 1 NaN NaN NaN ... 0 0 0 2734.687353 4258.980140 7772.212858 2 28 NaN NaN
13581 12568 2 NaN NaN NaN NaN 5 NaN 1 NaN ... 0 0 0 2734.687353 4258.980140 7772.212858 2 28 NaN 6.3750
13582 12568 3 NaN NaN NaN NaN 4 NaN NaN NaN ... 0 0 0 2734.687353 4258.980140 7772.212858 2 28 NaN NaN
13583 12569 1 NaN NaN NaN NaN 3 NaN NaN NaN ... 0 0 0 2580.967613 2925.167116 5075.164946 2 61 NaN NaN
13584 12569 2 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 2580.967613 2925.167116 5075.164946 2 61 NaN 6.3750
13585 12570 1 NaN NaN NaN NaN 3 NaN NaN NaN ... 0 0 0 5181.311509 6205.829154 11325.017623 2 40 NaN NaN
13586 12570 2 NaN NaN NaN NaN 3 NaN NaN NaN ... 0 0 0 5181.311509 6205.829154 11325.017623 2 40 NaN NaN
13587 12570 3 NaN NaN NaN NaN 3 NaN NaN NaN ... 0 0 0 5181.311509 6205.829154 11325.017623 2 40 NaN NaN
13588 12571 1 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 4670.540953 5795.692880 6269.200989 1 78 NaN 6.1875
13589 12571 2 NaN NaN NaN NaN 3 NaN NaN NaN ... 0 0 0 4670.540953 5795.692880 6269.200989 1 78 NaN NaN
13590 12571 3 NaN NaN NaN NaN 3 NaN NaN NaN ... 0 0 0 4670.540953 5795.692880 6269.200989 1 78 NaN NaN
13591 12571 4 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 4670.540953 5795.692880 6269.200989 1 78 NaN 7.5000
13592 12571 5 NaN NaN NaN NaN 6 NaN 1 NaN ... 0 0 0 4670.540953 5795.692880 6269.200989 1 78 NaN 7.5000

13593 rows × 244 columns

Select live births, then make a CDF of totalwgt_lb.


In [4]:
import pandas as pd
preg_live = preg[ preg.outcome == 1]
preg_live_totalwgt_lb = pd.DataFrame(preg_live.totalwgt_lb.dropna())
preg_live_totalwgt_lb = preg_live_totalwgt_lb.sort(['totalwgt_lb']).reset_index(drop=True)
cdf_dict = {}

for index, row in preg_live_totalwgt_lb.iterrows():
    cdf_dict[row[0]] = (index/ float(len(preg_live_totalwgt_lb) - 1))
    

cdf = pd.DataFrame(cdf_dict.keys(), columns = ['values'])
cdf['p_rank'] = cdf_dict.values()
cdf =  cdf.sort(['p_rank']).reset_index(drop=True)
cdf


Out[4]:
values p_rank
0 0.1250 0.000000
1 0.3125 0.000111
2 0.4375 0.000221
3 0.5625 0.000332
4 0.6250 0.000553
5 0.9375 0.000664
6 1.0000 0.000775
7 1.0625 0.001107
8 1.1250 0.001328
9 1.1875 0.001439
10 1.2500 0.001660
11 1.3125 0.001992
12 1.3750 0.002656
13 1.4375 0.002766
14 1.5000 0.003209
15 1.5625 0.003541
16 1.6250 0.004094
17 1.6875 0.004205
18 1.7500 0.004648
19 1.8125 0.004758
20 1.8750 0.005090
21 2.0000 0.005533
22 2.0625 0.005975
23 2.1250 0.006307
24 2.1875 0.007082
25 2.2500 0.007525
26 2.3125 0.007967
27 2.3750 0.008299
28 2.4375 0.008631
29 2.5000 0.009184
... ... ...
154 10.4375 0.990926
155 10.5000 0.992254
156 10.5625 0.992586
157 10.6250 0.993029
158 10.6875 0.993914
159 10.7500 0.994246
160 10.8125 0.994799
161 10.8750 0.995131
162 10.9375 0.995242
163 11.0000 0.996238
164 11.0625 0.996348
165 11.1250 0.996570
166 11.1875 0.996902
167 11.3750 0.997012
168 11.4375 0.997234
169 11.5000 0.997344
170 11.5625 0.997566
171 11.6250 0.997676
172 11.6875 0.997787
173 11.7500 0.998008
174 11.9375 0.998119
175 12.0000 0.998672
176 12.1875 0.998893
177 12.3750 0.999115
178 12.5000 0.999225
179 13.0000 0.999336
180 13.5000 0.999447
181 13.7500 0.999557
182 14.0000 0.999889
183 15.4375 1.000000

184 rows × 2 columns

Display the CDF.


In [5]:
cdf.plot(x = 'values', y = 'p_rank', legend = False)


Out[5]:
<matplotlib.axes._subplots.AxesSubplot at 0x10ae6a6d0>

In [6]:
# lambad function to create cdf of any series
create_cdf_pandas = lambda data: data.value_counts().sort_index().cumsum()*1./len(data)

create_cdf_pandas(preg_live.totalwgt_lb.dropna()).plot(legend=True, label = 'cdf')

ser = create_cdf_pandas(preg_live.totalwgt_lb.dropna())


Find out how much you weighed at birth, if you can, and compute CDF(x).


In [7]:
cdf_value = lambda cdf, num : cdf[num] if num in cdf else cdf[min(cdf.keys(), key=lambda k: abs(k-num))]

cdf_value(ser.to_dict(), 7.5)


Out[7]:
0.57756140739101569

If you are a first child, look up your birthweight in the CDF of first children; otherwise use the CDF of other children.

Compute the percentile rank of your birthweight


In [192]:
#first_live = preg[preg.pregordr == 1 & preg.outcome == 1]
#first_live
df = preg.query('pregordr == 1 and outcome == 1')
ser = create_cdf_pandas(df.totalwgt_lb.dropna())
cdf_value(ser.to_dict(), 7.5)


Out[192]:
0.60216671682214862

Compute the median birth weight by looking up the value associated with p=0.5.


In [8]:
#find_median_value = lambda cdf_dict: [key for key, value in cdf_dict.iteritems() if value == 0.5 
#find_median_value(ser.to_dict())
cdf_dict = ser.to_dict()
def get_key(cdf_dict, p):
    if 0.5 in cdf_dict.values():
        for key, value in cdf_dict.iteritems():
            if value == p:
                return key
    else:
        return min(cdf_dict, key=lambda y:abs(float(cdf_dict[y]) - p))
print get_key(cdf_dict, 0.5)


7.3125

Compute the interquartile range (IQR) by computing percentiles corresponding to 25 and 75.


In [38]:
ser = create_cdf_pandas(preg_live.totalwgt_lb.dropna())
cdf_dict = ser.to_dict()

print get_key(cdf_dict, 0.25)
print get_key(cdf_dict, 0.75)


6.4375
8.0625

Make a random selection from cdf.


In [39]:
pd.DataFrame(ser).sample()


Out[39]:
0
7.5 0.577561

Draw a random sample from cdf.


In [40]:
pd.DataFrame(ser).sample(n = 10)


Out[40]:
0
14.0000 0.999889
10.6250 0.993029
10.9375 0.995242
6.7500 0.325957
1.3125 0.002102
8.6875 0.873977
5.3750 0.080438
1.6250 0.004204
9.0000 0.924873
3.7500 0.020137

Draw a random sample from cdf, then compute the percentile rank for each value, and plot the distribution of the percentile ranks.


In [53]:
df1 = pd.DataFrame(ser).sample(n = 1000, replace = True)
df1.reset_index(inplace=True)
df1.columns = ['val', 'prob']
df1_sample = df1.val
sample_ser = create_cdf_pandas(df1_sample)
#print len(ser)
#print df1_sample.value_counts()
#print df1_sample.value_counts().sort_index()
#print pd.DataFrame(df1, columns = ['val', 'prob'])
sample_ser.plot()


Out[53]:
<matplotlib.axes._subplots.AxesSubplot at 0x10d50cb10>

Generate 1000 random values using random.random() and plot their PMF.


In [65]:
import random
t = [random.random() for _ in range(1000)]
ser_pmf = pd.Series(t).value_counts().sort_index()/ len(t)
ser_pmf.plot()


Out[65]:
<matplotlib.axes._subplots.AxesSubplot at 0x10d596690>

Assuming that the PMF doesn't work very well, try plotting the CDF instead.


In [66]:
cdf_ser = create_cdf_pandas(pd.Series(t))

In [68]:
cdf_ser.plot()


Out[68]:
<matplotlib.axes._subplots.AxesSubplot at 0x10d51f350>

In [ ]: