In [18]:

    

from simpleabc import simple_abc
import simple_model
import numpy as np
import pickle 
import pylab as plt
from scipy import stats
import time
%matplotlib inline
plt.style.use('ggplot')


    

In [30]:

    

reload(simple_model)


    

    
Out[30]:





<module 'simple_model' from 'simple_model.py'>

In [31]:

    

np.random.seed(914)

steps = 3
eps = 0.25
min_part = 10

#stars = pickle.load(file('stars.pkl'))
stars = pickle.load(file('stars_trimmed.pkl'))
#obs = pickle.load(file('data.pkl'))

model = simple_model.MyModel(stars)

#theta = (mututal inclination, eccentricity, planet number)

model.set_prior([stats.uniform(0, 90.0),
                stats.uniform(0,20)])



theta = (2.0,10)
obs = model.generate_data(theta)
model.set_data(obs)





n_procs = [1, 2, 3, 4, 5, 6, 7, 8]

start = time.time()
OT = simple_abc.pmc_abc(model, obs, epsilon_0=eps, min_particles=min_part, steps=steps,
                        target_epsilon=eps, parallel=False)
end = time.time()
print 'Serial took {}s'.format(end - start)
out_pickle = file('simptest.pkl', 'w')
pickle.dump(OT, out_pickle)
out_pickle.close()


    

    




<type 'tuple'>
0 0.25
1 0.0903982622764
Effective sample size(s): [[ 11.  11.]]
2 0.0790279005374
Effective sample size(s): [[ 11.   7.]]
Serial took 102.543305874s

In [32]:

    

for P in OT:
    P[0] = P[0].T
    plt.figure(figsize=(15,3))
    plt.subplot(121)
    ker = stats.gaussian_kde([x[0] for x in P[0]])
    plt.hist([x[0] for x in P[0]],normed=True, alpha=0.2)
    x = np.linspace(0.0,90,1000)
    plt.plot(x,ker(x))
    plt.axvline(theta[0])
    plt.xlabel(r"$\sigma_{mi}$")
    plt.subplot(122)
    ker = stats.gaussian_kde([x[1] for x in P[0]])
    plt.hist([x[1] for x in P[0]],normed=True, alpha=0.2)
    x = np.linspace(0,20,1000)
    plt.plot(x,ker(x))
    plt.axvline(theta[1])
    plt.xlabel(r"$\lambda$")


    

In [2]:

    

DT = pickle.load(file('demo.pkl'))


    

In [3]:

    

thetad = (2.0, 0.05, 10)
for P in DT:
    plt.figure(figsize=(15,3))
    plt.subplot(131)
    ker = stats.gaussian_kde([x[0] for x in P[0]])
    plt.hist([x[0] for x in P[0]],normed=True, alpha=0.2)
    x = np.linspace(0.0,90,1000)
    plt.plot(x,ker(x))
    plt.axvline(thetad[0])
    plt.xlabel(r"$\sigma_{mi}$")
    plt.subplot(132)
    ker = stats.gaussian_kde([x[1] for x in P[0]])
    plt.hist([x[1] for x in P[0]],normed=True, alpha=0.2)
    x = np.linspace(0.0,1.0,1000)
    plt.plot(x,ker(x))
    plt.axvline(thetad[1])
    plt.xlabel(r"$\sigma_{e}$")
    plt.subplot(133)
    ker = stats.gaussian_kde([x[2] for x in P[0]])
    plt.hist([x[2] for x in P[0]],normed=True, alpha=0.2)
    x = np.linspace(0,35,1000)
    plt.plot(x,ker(x))
    plt.axvline(thetad[2])
    plt.xlabel(r"$\lambda$")


    

In [5]:

    

epsilon = [DT[i][4] for i in range(0,20)]
plt.plot(range(0,20),epsilon, 'o-')
accept = [DT[i][2]/float(DT[i][3]) for i in range(0,20)]
plt.plot(range(0,20),accept, 'o-')


    

    
Out[5]:





[<matplotlib.lines.Line2D at 0x10f72d910>]






    

In [13]:

    




    

    
Out[13]:





([[39.81183953660224, 0.17249596105402232, 14.61474273766613], [39.6598740454627, 0.18772377679374885, 10.016040614078442], [47.483679253949674, 0.18334077497836965, 7.864872626610342], [78.31508908045609, 0.03199776353614603, 9.36693274906407], [83.365683518924, 0.11020693940364323, 17.008628237965787], [11.285233349332128, 0.26787275344339523, 14.358438776632887], [22.8802079826819, 0.1869143076490163, 17.178789955116407], [47.231436130154556, 0.1882786288254532, 6.700995187095957], [39.43517466492983, 0.2676446498218684, 19.390263413399822], [48.628469834246644, 0.313236615412126, 15.808245170373008], [29.877438489361207, 0.3107168821809947, 6.27956095718039], [38.886712951081435, 0.33235183238971155, 15.323385713107477], [81.2959354088329, 0.1553327623075893, 19.86082643267811], [19.706581919177015, 0.1941846628711823, 9.992025496564397], [84.8796752204316, 0.12497365746940303, 11.77710022377886], [73.0198351379791, 0.10000703775606257, 18.969379774716508], [53.371373371323145, 0.18391020947384895, 9.392254325186393], [77.04037241462844, 0.0904151269086566, 12.205161427359384], [62.44111470678974, 0.09969988515854578, 8.130496786865528], [50.73842512076609, 0.2458901149198809, 10.814827717254506], [76.7975326970149, 0.0787250778756724, 9.652137637300473], [15.177368809219242, 0.001497379243232122, 17.027328079757385], [3.8682072909205987, 0.3397157350687834, 7.597304891256744], [34.497762940662064, 0.21542389015284957, 9.519386415340108], [31.99692648602532, 0.1578835647463298, 13.540119088159168], [88.42473678643398, 0.14050928822071007, 11.052116902964306], [57.96667059119036, 0.2591336809966399, 9.256763090790303], [33.78865980828721, 0.1320095500601326, 15.306323762333614], [37.4219782625993, 0.1726569241418341, 19.250689142843292], [41.750033027841745, 0.3113703023788692, 13.392296443197527], [45.603649100339645, 0.3067346265246419, 14.211783689149826], [66.62600958334183, 0.0035167735845559323, 9.011838785008448], [66.210016443036, 0.02593037465487591, 11.167416324811594], [24.43195966498331, 0.2318877230422346, 12.026894754402566], [27.944963791633743, 0.1583342877537245, 13.224326419420343], [2.2438308681027754, 0.060660732177442944, 10.04422731048371], [67.69047349878016, 0.013021936774152199, 14.502588990676117], [62.75462016923507, 0.18821225257451613, 16.79280941360221], [17.166374703160987, 0.13097478495449288, 12.279177517289341], [54.72264984448811, 0.2838802528081701, 8.352960460000467], [62.33252190423285, 0.09221979439399097, 18.44283585662027], [85.12591235092576, 0.25115036116675016, 18.843051347491006], [20.05548164027834, 0.20799608264276903, 18.654452141803198], [13.11203809584254, 0.06254782917411761, 5.709753298870801], [52.691877832549324, 0.30308356528192915, 18.2818808622972], [83.88005710053177, 0.17000401858985648, 8.055859406629054], [49.74281043548792, 0.23473865572804464, 19.26501664151681], [17.933021916195194, 0.09783053454136925, 7.843298500555194], [53.466425342031556, 0.032347512423414826, 17.256724315709814], [2.4881937213791705, 0.09225684213256957, 2.671261542615624], [45.25348463246463, 0.075438937143389, 13.453964541760868], [5.787251869238358, 0.005716623403865806, 12.714206119628733], [20.95254391038876, 0.1361429069716733, 16.579996890844154], [66.51306783080503, 0.33922711799922944, 16.636840184851906], [48.57358040199928, 0.015445583159260923, 10.512358867839724], [4.550870364474604, 0.3063787350217134, 10.651254419182749], [39.059298049131634, 0.2288931167549657, 8.922888653270377], [79.63032344265048, 0.328212962563291, 19.735000274118647], [8.25395306238835, 0.12785751684356406, 13.016260065517045], [47.320211662444564, 0.08146851063937255, 18.119929376489697], [48.56182931236124, 0.3012937760804558, 12.934869063926717], [39.88772877831791, 0.3302834580908858, 18.465609509495426], [27.133318393181167, 0.20910260299731942, 10.042635070627009], [0.6518388296744593, 0.10763518240753644, 8.195462768438391], [87.09332308410194, 0.1889644803331263, 19.207410022783172], [69.66704431775506, 0.289415520250912, 15.137252598785654], [24.83208213417861, 0.23195565710611943, 6.423077359861078], [53.87825286449115, 0.2155052852168592, 13.043224294030582], [65.96013615461192, 0.28535289653123186, 11.348263751122944], [82.19490060642104, 0.31306357575806676, 18.5406098326644], [67.27615625871653, 0.1826792093094024, 11.498910957355815], [6.970577549444748, 0.19196081911359197, 10.81066768227528], [28.7512714331863, 0.007060828681669351, 6.797881467628601], [75.0550848666146, 0.017800049104868876, 11.874881203856258], [16.73068306066471, 0.007443559940987443, 5.728179828478021], [12.942716690880227, 0.13631866502432377, 14.248398517385167], [3.5458134698598953, 0.1018089660650422, 6.52237198347744], [17.954789704766323, 0.14571207913275075, 7.713173435669285], [16.593137175653208, 0.2503909244067938, 17.19791733250483], [89.55734033896046, 0.07994470503326412, 10.909596951209679], [65.01072332717288, 0.2594459099364874, 16.707503990142005], [32.74169886491271, 0.24081988528742637, 17.957504159559146], [56.9621101404683, 0.2212539108896141, 15.031952688067216], [43.70557672114724, 0.24125521224644264, 16.8185991213999], [33.7768207311514, 0.15454827919490732, 11.765442354664], [37.86158080300729, 0.16807144543578367, 19.822695273213363], [78.65888099971431, 0.13015073040122482, 8.795122675580732], [85.83752808912094, 0.18927953747634274, 14.07493678214675], [37.99550627080449, 0.006585270723693859, 11.044970798635191], [84.91014320725645, 0.14804094126290301, 19.365175123249244], [48.32251141146499, 0.12251434780862847, 15.005820132264788], [78.03850551903925, 0.3464200684903572, 12.566978257901983], [84.89295051928285, 0.2551350088397333, 10.243142368392498], [43.41062479239401, 0.26398553460036167, 11.700200942993186], [55.911890646990535, 0.055727302936845136, 14.213121855441193], [59.22589684509844, 0.24968710639098068, 13.642925544157789], [37.17701691037164, 0.24685445966146047, 17.228231748334434], [57.52252340371407, 0.19203297716436252, 19.977233108935554], [57.13667807807183, 0.20239290078014804, 7.778300877380467], [10.875156082456002, 0.002657671677603002, 15.614158747373772], [78.39779394953354, 0.043377945272894625, 17.510391991399366]], [0.085344245744205169, 0.14985244030144412, 0.20916827380669489, 0.20442043631745602, 0.10113529046272325, 0.21977832452002299, 0.14393200442551884, 0.23955516052608752, 0.16460144025679924, 0.19088426434103728, 0.2489367072708723, 0.21831034103179472, 0.087602066052472893, 0.088950367474356043, 0.16617847042297029, 0.081750479898558506, 0.19100045267919191, 0.15571186410641599, 0.23220807270361268, 0.16801768170792919, 0.20840298759908243, 0.19054014939944744, 0.21600914965221152, 0.14361704913566439, 0.074764480903188157, 0.17642225177457074, 0.20729975809545312, 0.071474096953652699, 0.099170492238067662, 0.18671461579394519, 0.18324971882941207, 0.21756002423034487, 0.17335683126150384, 0.10395661237100726, 0.07271625080940991, 0.198500542870508, 0.12050600286834039, 0.095397717769586207, 0.092606288267840536, 0.23327893001856745, 0.086484373044974708, 0.1291015917052403, 0.19229701647768754, 0.12971061807195344, 0.18642442456674327, 0.22525561363913996, 0.1187276757104475, 0.10733898955678271, 0.082757204845050739, 0.15870944757074343, 0.10471003458813177, 0.2276534096835377, 0.13214610034715005, 0.22431497373330428, 0.16324096474731664, 0.21812642986542033, 0.17268268425344088, 0.22290223840566278, 0.19563193423006064, 0.077365701156569863, 0.18916690487496576, 0.22715669648339987, 0.11321700817670598, 0.1440639795077642, 0.093120260693382598, 0.17101874446124476, 0.19128466625571597, 0.13035138159354445, 0.19593076832753967, 0.19582855970258672, 0.15562101802480643, 0.16570454780226077, 0.19740972614042085, 0.17303908449142025, 0.16499943745672277, 0.16163754693274046, 0.069961090648060636, 0.093758640610673322, 0.21237611070546378, 0.18069532199203991, 0.13691071278381145, 0.14175844706239624, 0.11992222167527991, 0.11933458547053251, 0.10033281764711811, 0.10682100441347442, 0.2144789847918758, 0.13538037193289601, 0.11708707124471578, 0.089342553374154296, 0.08777826194512961, 0.23704435401387489, 0.19775795260107462, 0.15645473731561815, 0.10920560940359059, 0.14804497257364463, 0.12974994856747771, 0.091263584089819896, 0.23112732810257772, 0.22112873592669088, 0.10340667122277983], 101.0, 428.0, 0.25)

In [ ]: