Adds best-fit Fmax distribution to noise plots and fix tables

Paper/AllSkyMC/generate_table.py

@@ -82,7 +82,7 @@ mcmc = pyfstat.MCMCFollowUpSearch(
     sftfilepath='{}/*{}*sft'.format(outdir, data_label),
     theta_prior=theta_prior,
     tref=tref, minStartTime=tstart, maxStartTime=tend,
-    nwalkers=nwalkers, ntemps=ntemps,
+    nwalkers=nwalkers, ntemps=ntemps, nsteps=[nsteps, nsteps]
     log10temperature_min=log10temperature_min)
-mcmc.run(Nsegs0=20, R=10)
-#mcmc.run(run_setup)
+#mcmc.run(Nsegs0=20, R=10)
+mcmc.run(run_setup)

Paper/AllSkyMCNoiseOnly/plot_data.py

@@ -5,6 +5,7 @@ import os
 from tqdm import tqdm
 from oct2py import octave
 import glob
+from scipy.stats import rv_continuous, chi2
 
 filenames = glob.glob("CollectedOutput/*.txt")
 
@@ -13,12 +14,13 @@ plt.style.use('paper')
 Tspan = 100 * 86400
 
 
-def maxtwoFinNoise(twoF, Ntrials):
-    F = twoF/2.0
-    alpha = (1 + F)*np.exp(-F)
-    a = Ntrials/2.0*F*np.exp(-F)
-    b = (1 - alpha)**(Ntrials-1)
-    return a*b
+class maxtwoFinNoise_gen(rv_continuous):
+    def _pdf(self, twoF, Ntrials):
+        F = twoF/2.0
+        alpha = (1 + F)*np.exp(-F)
+        a = Ntrials/2.0*F*np.exp(-F)
+        b = (1 - alpha)**(Ntrials-1)
+        return a*b
 
 
 df_list = []
@@ -31,6 +33,21 @@ df = pd.concat(df_list)
 print 'Number of samples = ', len(df)
 
 fig, ax = plt.subplots()
+
+maxtwoFinNoise = maxtwoFinNoise_gen(a=0)
+Ntrials_effective, loc, scale = maxtwoFinNoise.fit(df.twoF.values, floc=0, fscale=1)
+print 'Ntrials effective = {:1.2e}'.format(Ntrials_effective)
+twoFsmooth = np.linspace(0, df.twoF.max(), 1000)
+best_fit_pdf = maxtwoFinNoise.pdf(twoFsmooth, Ntrials_effective)
+ax.plot(twoFsmooth, best_fit_pdf, '-r')
+
+pval = 1e-6
+twoFsmooth_HD = np.linspace(
+    twoFsmooth[np.argmax(best_fit_pdf)], df.twoF.max(), 100000)
+best_fit_pdf_HD = maxtwoFinNoise.pdf(twoFsmooth_HD, Ntrials_effective)
+spacing = twoFsmooth_HD[1]-twoFsmooth_HD[0]
+print twoFsmooth_HD[np.argmin(np.abs(best_fit_pdf_HD - pval))], spacing
+
 ax.hist(df.twoF, bins=50, histtype='step', color='k', normed=True, linewidth=1)
 twoFsmooth = np.linspace(0, df.twoF.max(), 100)
 # ax.plot(twoFsmooth, maxtwoFinNoise(twoFsmooth, 8e5), '-r')

Paper/DirectedMC/generate_table.py

@@ -27,7 +27,7 @@ DeltaF1 = VF1 * np.sqrt(45/4.)/(np.pi*Tspan**2)
 
 depth = 100
 
-nsteps = 20
+nsteps = 25
 run_setup = [((nsteps, 0), 20, False),
              ((nsteps, 0), 7, False),
              ((nsteps, 0), 2, False),
@@ -71,6 +71,6 @@ mcmc = pyfstat.MCMCFollowUpSearch(
     sftfilepath='{}/*{}*sft'.format(outdir, data_label),
     theta_prior=theta_prior,
     tref=tref, minStartTime=tstart, maxStartTime=tend,
-    nwalkers=nwalkers, ntemps=ntemps,
+    nwalkers=nwalkers, ntemps=ntemps, nsteps=[nsteps, nsteps],
     log10temperature_min=log10temperature_min)
 mcmc.run(run_setup)

Paper/DirectedMCNoiseOnly/plot_data.py

@@ -5,6 +5,7 @@ import os
 from tqdm import tqdm
 from oct2py import octave
 import glob
+from scipy.stats import rv_continuous, chi2
 
 filenames = glob.glob("CollectedOutput/*.txt")
 
@@ -13,12 +14,13 @@ plt.style.use('paper')
 Tspan = 100 * 86400
 
 
-def maxtwoFinNoise(twoF, Ntrials):
-    F = twoF/2.0
-    alpha = (1 + F)*np.exp(-F)
-    a = Ntrials/2.0*F*np.exp(-F)
-    b = (1 - alpha)**(Ntrials-1)
-    return a*b
+class maxtwoFinNoise_gen(rv_continuous):
+    def _pdf(self, twoF, Ntrials):
+        F = twoF/2.0
+        alpha = (1 + F)*np.exp(-F)
+        a = Ntrials/2.0*F*np.exp(-F)
+        b = (1 - alpha)**(Ntrials-1)
+        return a*b
 
 
 df_list = []
@@ -32,8 +34,21 @@ print 'Number of samples = ', len(df)
 
 fig, ax = plt.subplots()
 ax.hist(df.twoF, bins=50, histtype='step', color='k', normed=True, linewidth=1)
-twoFsmooth = np.linspace(0, df.twoF.max(), 100)
-# ax.plot(twoFsmooth, maxtwoFinNoise(twoFsmooth, 2e3), '-r')
+
+maxtwoFinNoise = maxtwoFinNoise_gen(a=0)
+Ntrials_effective, loc, scale = maxtwoFinNoise.fit(df.twoF.values, floc=0, fscale=1)
+print 'Ntrials effective = {:1.2e}'.format(Ntrials_effective)
+twoFsmooth = np.linspace(0, df.twoF.max(), 1000)
+best_fit_pdf = maxtwoFinNoise.pdf(twoFsmooth, Ntrials_effective)
+ax.plot(twoFsmooth, best_fit_pdf, '-r')
+
+pval = 1e-6
+twoFsmooth_HD = np.linspace(
+    twoFsmooth[np.argmax(best_fit_pdf)], df.twoF.max(), 100000)
+best_fit_pdf_HD = maxtwoFinNoise.pdf(twoFsmooth_HD, Ntrials_effective)
+spacing = twoFsmooth_HD[1]-twoFsmooth_HD[0]
+print twoFsmooth_HD[np.argmin(np.abs(best_fit_pdf_HD - pval))], spacing
+
 ax.set_xlabel('$\widetilde{2\mathcal{F}}$')
 ax.set_xlim(0, 60)
 fig.tight_layout()

Paper/allsky_setup_run_setup.tex

 \begin{tabular}{c|cccccc}
 Stage & $\Nseg$ & $\Tcoh^{\rm days}$ &$\Nsteps$ & $\V$ & $\Vsky$ & $\Vpe$ \\ \hline
-0 & 20 & 5.0 & 100 & $4{\times}10^{2}$ & 8.0 & 60.0 \\
-1 & 11 & 9.1 & 100 & $4{\times}10^{3}$ & 20.0 & $2{\times}10^{2}$ \\
-2 & 6 & 16.7 & 100 & $4{\times}10^{4}$ & 70.0 & $6{\times}10^{2}$ \\
-3 & 3 & 33.3 & 100 & $3{\times}10^{5}$ & $1{\times}10^{2}$ & $2{\times}10^{3}$ \\
-4 & 1 & 100.0 & 100,100 & $2{\times}10^{6}$ & $2{\times}10^{2}$ & $1{\times}10^{4}$ \\
+0 & 20 & 5.0 & 50 & $4{\times}10^{2}$ & 8.0 & 60.0 \\
+1 & 11 & 9.1 & 50 & $4{\times}10^{3}$ & 20.0 & $2{\times}10^{2}$ \\
+2 & 6 & 16.7 & 50 & $4{\times}10^{4}$ & 70.0 & $6{\times}10^{2}$ \\
+3 & 3 & 33.3 & 50 & $3{\times}10^{5}$ & $1{\times}10^{2}$ & $2{\times}10^{3}$ \\
+4 & 1 & 100.0 & 50,50 & $2{\times}10^{6}$ & $2{\times}10^{2}$ & $1{\times}10^{4}$ \\
 \end{tabular}

Paper/directed_setup_run_setup.tex

 \begin{tabular}{c|cccc}
 Stage & $\Nseg$ & $\Tcoh^{\rm days}$ &$\Nsteps$ & $\Vpe$ \\ \hline
-0 & 20 & 5.0 & 20 & 10.0 \\
-1 & 7 & 14.3 & 20 & $1{\times}10^{2}$ \\
-2 & 2 & 50.0 & 20 & $1{\times}10^{3}$ \\
-3 & 1 & 100.0 & 20,20 & $2{\times}10^{3}$ \\
+0 & 20 & 5.0 & 25 & 10.0 \\
+1 & 7 & 14.3 & 25 & $1{\times}10^{2}$ \\
+2 & 2 & 50.0 & 25 & $1{\times}10^{3}$ \\
+3 & 1 & 100.0 & 25,25 & $2{\times}10^{3}$ \\
 \end{tabular}

Paper/macros.tex

-\def\AllSkyMCNoiseN{10000}
+\def\AllSkyMCNoiseN{1.0{\times}10^{4}}
 \def\AllSkyMCNoiseOnlyMaximum{59.8}
 \def\BasicExampleDeltaFone{1.0{\times}10^{-13}}
 \def\BasicExampleDeltaFzero{1.0{\times}10^{-7}}
@@ -12,5 +12,5 @@
 \def\BasicExamplehzero{1.0{\times}10^{-24}}
 \def\BasicExamplenburn{50.0}
 \def\BasicExamplenprod{50.0}
-\def\DirectedMCNoiseN{10000}
+\def\DirectedMCNoiseN{1.0{\times}10^{4}}
 \def\DirectedMCNoiseOnlyMaximum{52.4}