diff --git a/src/seismic_hazard_forecasting.py b/src/seismic_hazard_forecasting.py
index 429a7c2..e770159 100644
--- a/src/seismic_hazard_forecasting.py
+++ b/src/seismic_hazard_forecasting.py
@@ -100,6 +100,7 @@ def plot_results(act_rate, bin_edges, bin_edges_dt, rt, boundaries,
     ax.patch.set_visible(False)
     
     fig.tight_layout()
+    plt.savefig("activity_rate.png", dpi=600)
     plt.show()
 
 def bootstrap_forecast(data):
@@ -613,6 +614,10 @@ verbose: {verbose}")
         sorted_pairs = sorted(zip(datenum_data, mag_data), key=lambda x: x[0])
         datenum_data, mag_data = map(list, zip(*sorted_pairs))
 
+        #-------split the data into bins and apply BEAST for changepoint detection--------------------
+        act_rate, bin_counts, bin_edges, bin_edges_dt, out, cps, rt, boundaries, bin_dur, time_unit = bins_and_beast(
+            np.array(datenum_data), time_unit, time_range, multiplicator)
+        
         #------Forecasted rate is taken from BEAST or is equal to last value if no changepoints detected-----
         if len(cps) > 0:
             rate_forecast = rt[-1]