{"cells":[{"metadata":{"_uuid":"fee8782af956b37d4a43de6736bd2b8eaa7f07b9"},"cell_type":"markdown","source":"Original code was adopted from https://github.com/aerdem4/kaggle-plasticc/blob/master/bazin.py . Thanks AhmetErdem for sharing the code!\nThis is a notebook to demonstrate how some simple optimisation tricks can help reduce time spent on feature enginnering and speed up your model improvement iterations.\n\n### original code:"},{"metadata":{"trusted":true,"_uuid":"d2ab40d4bf7f27394bda1b680f3744d2c56dfaf6"},"cell_type":"code","source":"#### import argparse\nimport multiprocessing\nimport numpy as np\nimport pandas as pd\nfrom scipy.optimize import least_squares\nimport time\nimport warnings\nimport numba\nimport dask.dataframe as dd\nimport pdb\nfrom dask.distributed import LocalCluster, Client\n\nNUM_PARTITIONS = 100\nLOW_PASSBAND_LIMIT = 3\nFEATURES = [\"A\", \"B\", \"t0\", \"tfall\", \"trise\", \"cc\", \"fit_error\", \"status\", \"t0_shift\"]","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"# bazin, errorfunc and fit_scipy are developed using:\n# https://github.com/COINtoolbox/ActSNClass/blob/master/examples/1_fit_LC/fit_lc_parametric.py\ndef bazin(time, low_passband, A, B, t0, tfall, trise, cc):\n    X = np.exp(-(time - t0) / tfall) / (1 + np.exp((time - t0) / trise))\n    return (A * X + B) * (1 - cc * low_passband)\n\n\ndef errfunc(params, time, low_passband, flux, weights):\n    return abs(flux - bazin(time, low_passband, *params)) * weights\n\n\ndef fit_scipy(time, low_passband, flux, flux_err):\n    time -= time[0]\n    sn = np.power(flux / flux_err, 2)\n    start_point = (sn * flux).argmax()\n\n    t0_init = time[start_point] - time[0]\n    amp_init = flux[start_point]\n    weights = 1 / (1 + flux_err)\n    weights = weights / weights.sum()\n    guess = [0, amp_init, t0_init, 40, -5, 0.5]\n\n    result = least_squares(errfunc, guess, args=(time, low_passband, flux, weights), method='lm')\n    result.t_shift = t0_init - result.x[2]\n\n    return result\n\n\ndef yield_data(meta_df, lc_df):\n    cols = [\"object_id\", \"mjd\", \"flux\", \"flux_err\", \"low_passband\"]\n    for i in range(NUM_PARTITIONS):\n        yield meta_df[(meta_df[\"object_id\"] % NUM_PARTITIONS) == i][\"object_id\"].values, \\\n              lc_df[(lc_df[\"object_id\"] % NUM_PARTITIONS) == i][cols]\n\n\ndef get_params(object_id_list, lc_df, result_queue):\n    results = {}\n    for object_id in object_id_list:\n        light_df = lc_df[lc_df[\"object_id\"] == object_id]\n        try:\n            result = fit_scipy(light_df[\"mjd\"].values, light_df[\"low_passband\"].values,\n                               light_df[\"flux\"].values, light_df[\"flux_err\"].values)\n            results[object_id] = np.append(result.x, [result.cost, result.status, result.t_shift])\n        except Exception as e:\n            print(e)\n            results[object_id] = None\n    result_queue.put(results)\n\n\ndef parallelize(meta_df, df):\n    pool_size = multiprocessing.cpu_count() * 2\n    pool = multiprocessing.Pool(processes=pool_size)\n\n    manager = multiprocessing.Manager()\n    result_queue = manager.Queue()\n\n    for m, d in yield_data(meta_df, df):\n        pool.apply_async(get_params, (m, d, result_queue))\n\n    pool.close()\n    pool.join()\n\n    return [result_queue.get() for _ in range(NUM_PARTITIONS)]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"e79a090ae2abeae177377a5e6c6d201fe80ac873"},"cell_type":"code","source":"meta_df = pd.read_csv('../input/training_set_metadata.csv')\nlc_df = pd.read_csv('../input/training_set.csv')\nlc_df[\"low_passband\"] = (lc_df[\"passband\"] < LOW_PASSBAND_LIMIT).astype(int)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"3ee5ef55427fc7a5a042701bbd174a55b51bc097","scrolled":true},"cell_type":"code","source":"warnings.simplefilter('ignore', RuntimeWarning)\nstart = time.time()\nresult_list = parallelize(meta_df, lc_df)\ncalc_fin = time.time()\nprint(f'calculation time: {calc_fin - start}')","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"84f02387bdcf6584898b350e68b4b3b27096e357"},"cell_type":"markdown","source":"### Numba.jit trick\nUsing a simple @numba.jit with minor code modifications can sometimes drastically speed up code execution speed, especially if you have repeated function calls or loops. Here it is able to reduce run time to less than 50%."},{"metadata":{"trusted":true,"_uuid":"49c2c87729f2471092555b8aa441ef77ee14cf74"},"cell_type":"code","source":"# optimise 1 - using numba\n@numba.jit(nopython=True)  # added\ndef bazin(time, low_passband, A, B, t0, tfall, trise, cc):\n    X = np.exp(-(time - t0) / tfall) / (1. + np.exp((time - t0) / trise))\n    return (A * X + B) * (1. - cc * low_passband)\n\n@numba.jit(nopython=True)  # added\ndef errfunc(params, time, low_passband, flux, weights):\n    A, B, t0, tfall, trise, cc = params\n    return np.abs(flux - bazin(time, low_passband, A, B, t0, tfall, trise, cc)) * weights","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"900eb393d8b4a02b84e587d20965c4050eb9975f"},"cell_type":"code","source":"warnings.simplefilter('ignore', RuntimeWarning)\nstart = time.time()\nresult_list = parallelize(meta_df, lc_df)\ncalc_fin = time.time()\nprint(f'calculation time: {calc_fin - start}')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"d7d97aa10291132b8d124e0643eb845fd8f42426","scrolled":true},"cell_type":"code","source":"start = time.time()\nfinal_result = {}\nfor res in result_list:\n    final_result.update(res)\n\nfor index, col in enumerate(FEATURES):\n    meta_df[col] = meta_df[\"object_id\"].apply(lambda x: final_result[x][index])\nend = time.time()\nprint(f'feature collection time: {end - start}')","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"0318724b479097e9970552f12385c17be93e5197"},"cell_type":"markdown","source":"### Using Dask\nUsing dask does not really speed up code execution (it might even slow your code if your multiprocessing is optimised), but it saved you a lot of effort trying to reshape, split and combine yoru data."},{"metadata":{"trusted":true,"_uuid":"2bcfadc9c83afd279c613a64ec67c5885a4ecabb"},"cell_type":"code","source":"# optimise 2 - using numba, dask\ndef collect_results(chunk):\n    result = fit_scipy(time=chunk['mjd'].values,\n                       low_passband=chunk['low_passband'].values, \n                       flux=chunk['flux'].values, \n                       flux_err=chunk['flux_err'].values)\n    return pd.Series(np.append(result.x, [result.cost, result.status, result.t_shift]), index=FEATURES)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"15af01afe1a779a04cc2e93f2696bf780a356535"},"cell_type":"code","source":"cluster = LocalCluster(n_workers=4, processes=True, scheduler_port=0,\n                           diagnostics_port=8787)\nclient = Client(cluster)\nlc_df_dask = dd.from_pandas(lc_df, npartitions=1).set_index('object_id').repartition(npartitions=100)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"scrolled":false,"_uuid":"31072f8ff43997e5fb16498bd83060c007a1a082"},"cell_type":"code","source":"warnings.simplefilter('ignore', RuntimeWarning)\nstart = time.time()\n# result_list = lc_df.head(2000).groupby('object_id').apply(collect_results)\nresult_list = lc_df_dask.groupby('object_id').apply(\n    collect_results,\n    meta=dict(A='f8', B='f8', t0='f8', tfall='f8', trise='f8', \n             cc='f8', fit_error='f8', status='f8', t0_shift='f8')).compute()\nend = time.time()\nprint(end - start)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"51b35d39ce4cbc0b4aa45e0ab695f2281368257f"},"cell_type":"code","source":"result_list.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"0c71d815b4822e7b46176a14e60a33512ef2b7b4"},"cell_type":"markdown","source":"As we can see above, the code using dask is much more readable while not that much slower than manual multiprocessing implementation. The result is already in nice Pandas dataframe format."}],"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"name":"python","version":"3.6.6","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"}},"nbformat":4,"nbformat_minor":1}