{"cells":[{"metadata":{},"cell_type":"markdown","source":"### 1. test plasticc data structure"},{"metadata":{},"cell_type":"markdown","source":"the question is how to define the length of the time sequence while feeding the training data ?\n\nTake 'simple' format as an example. If we need to consider time dilation? Assume that we only use the luminosity variability data, without the knowledge of redshift. We have to set the identical length of our train and test light curves. For quasar classification, the point is not about where is the start, but how to extract the most useful part of their original data.\n\nGiven the plasticc data structure, for each mjd, about one band's data is given and other bands' are blank. Then how to consider these vacancies?"},{"metadata":{"_uuid":"35aada01-d9a4-4ea0-b0da-6f667e34b2cb","_cell_guid":"cdc983dd-339d-4959-9216-7dd106b3031a","trusted":true,"scrolled":false},"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nimport matplotlib.pyplot as plt \n\ndef diff_mjd(x):\n    return x.max()-x.min()\n\ndef early_mjd(x):\n    return x.min()\n\ndef late_mjd(x):\n    return x.max()\n\ndata = pd.read_csv('../input/plasticc-converted-datasets/training_set_converted.csv')\ndata = data[['mjd','id']]\ncombine = data.groupby('id')\n\ndiff = combine.agg(diff_mjd)\nearly = combine.agg(early_mjd)\nlate = combine.agg(late_mjd)\n\nmjds = data['mjd']\nprint('latest mjd: ',max(mjds))\nprint('earlies mjd: ',min(mjds))\n\nn_bins = 5\nplt.figure(figsize=(12, 8), dpi=100)\nplt.hist(diff['mjd'],bins=n_bins)\nplt.xlabel('mjd difference', fontsize=16)\nplt.ylabel('number',fontsize=16)\nplt.show()\n\nplt.figure(figsize=(12, 8), dpi=100)\nplt.hist(early['mjd'],bins=n_bins)\nplt.xlabel('mjd early', fontsize=16)\nplt.ylabel('number',fontsize=16)\nplt.show()\n\nplt.figure(figsize=(12, 8), dpi=100)\nplt.hist(late['mjd'],bins=n_bins)\nplt.xlabel('mjd late', fontsize=16)\nplt.ylabel('number',fontsize=16)\nplt.show()\nprint('max length of the light curve: ',max(diff['mjd']), '\\n min length: ', min(diff['mjd']))\nprint('training objects number: ',len(combine))\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### 2. light curve analysis"},{"metadata":{"trusted":true},"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nimport matplotlib.pyplot as plt \nimport random\n\n\ndata = pd.read_csv('../input/plasticc-converted-datasets/test_set_converted_batch2_2.csv')\ncombine = data.groupby('id')\nid_list = []\nfor _id,_group in combine:\n    id_list.append(_id)\n\n# random select one object's light curve\n# ram_id = random.choice(id_list)\nram_id =2744261\nprint(ram_id)\nobj_data = data[(data['id']==ram_id)]\n\n\ntime = obj_data['mjd']\n\n\nplt.figure(figsize = (20,10))\n\nplt.errorbar(time, obj_data['g'],yerr=obj_data['g_err'], alpha = 0.6,fmt='o',ecolor='green',color='green',elinewidth=2,capsize=3,ms=5,label = 'g')\nplt.errorbar(time, obj_data['u'],yerr=obj_data['u_err'], alpha = 0.6,fmt='o',ecolor='blue',color='blue',elinewidth=2,capsize=3,ms=5,label = 'u')\nplt.errorbar(time, obj_data['r'],yerr=obj_data['r_err'], alpha = 0.6,fmt='o',ecolor='red',color='red',elinewidth=2,capsize=3,ms=5,label = 'r')\nplt.errorbar(time, obj_data['i'],yerr=obj_data['i_err'], alpha = 0.6,fmt='o',ecolor='orange',color='orange',elinewidth=2,capsize=3,ms=5,label = 'i')\nplt.errorbar(time, obj_data['z'],yerr=obj_data['z_err'], alpha = 0.6,fmt='o',ecolor='black',color='black',elinewidth=2,capsize=3,ms=5,label = 'z')\nplt.errorbar(time, obj_data['y'],yerr=obj_data['y_err'], alpha = 0.6,fmt='o',ecolor='purple',color='purple',elinewidth=2,capsize=3,ms=5,label = 'y')\n\nplt.legend(loc=\"upper left\",fontsize=24)\nplt.xlabel(\"Time (days)\",fontsize=24)\nplt.ylabel(\"Flux\",fontsize=24)\nplt.xticks(fontsize=20)\nplt.yticks(fontsize=20)\nplt.show()\n# plt.gca().invert_yaxis()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"obj_data","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"It is easy to recognize that all light curves have two-three big gaps, which is caused by the switch of telescope's observe erea. From these, we could still build 'group' format, but with more flexible boundary dates.\nThe next step is to calculate the suitable gap. Btw, the boundary limitation should also be improved in our original codes. "},{"metadata":{"trusted":true},"cell_type":"code","source":"# compare with the original data\ndata = pd.read_csv('../input/PLAsTiCC-2018/test_set_batch2.csv')\nram_id =2744261\nobj_data = data[(data['object_id']==ram_id)]\nobj_data\n# okay,my convert file has problem! I remember in STATA, there is a method to convert the values of a column to header names. \n# That would be the quickest way to convert files without errors!!!","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### 3. mjd gaps analysis"},{"metadata":{"trusted":true},"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nimport matplotlib.pyplot as plt \nimport random\n\n\ndata = pd.read_csv('../input/plasticc-converted-datasets/training_set_converted.csv')\ncombine = data.groupby('id')\nid_list = []\nfor _id,_group in combine:\n    id_list.append(_id)\n\ndelta_time = []\nobservation_time = []\nfor i in id_list:\n    ind_list = data[data.id == i].index.tolist()\n    mjd_exit_list = data['mjd'][ind_list].tolist()\n    m = 0\n    while m<len(mjd_exit_list)-1:\n        delta_time.append(mjd_exit_list[m+1]-mjd_exit_list[m])\n        m +=1\n    observation_time.append(mjd_exit_list[-1]-mjd_exit_list[0])\n\n\n# plot\n\nplt.figure(figsize = (8,8))\n# plt.subplot(1,2,1)\nplt.hist(x=delta_time,bins=10**np.arange(0,3,0.1))# bins=10**np.arange(0,3,0.1)\nplt.xlabel(r'$\\Delta$ t (days)', fontsize=20)\nplt.ylabel('Number', fontsize=20)\nplt.xscale('log')\nplt.savefig('delta_time.png')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"from the plot we could see, delta t larger than 10^2 could be the gap between each season. In this way, we could estimate the suitable gap for each light curve."},{"metadata":{"trusted":true},"cell_type":"code","source":"# plot the observation times for each object\nplt.figure(figsize = (8,8))\nobj_list = combine['id'].count().tolist()\nplt.hist(x=obj_list)\nplt.xlabel('data length', fontsize=20)\nplt.ylabel('Number', fontsize=20)\nplt.savefig('data_length.png')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### 4. find suitable gap length"},{"metadata":{"trusted":true},"cell_type":"code","source":"def remove_alone_mjd(id_list = list, data = list, check_delta = 300, min_size = 3):\n\t'''\n\tThis function is used to remove the data points whose mjd is far away from other data points.\n\tThese alone data pionts is not useful for group format and season format input data.\n\tInput:\n\t- id_list: the id list of all objects\n\t- data: the data from the preprocessed file\n\t- check_delta: if the difference between two neighboring data points' mjd is larger than this value, the earlier one will be removed.\n\t- min_size: the minimal size of the group before padding\n\tReturns:\n\t- data: the modified data\n\t'''\n\n\trecord_row = []\n\n\tfor i in id_list:\n\t\tind_list = data[data.object_id == i].index.tolist()\n\t\tmjd_exit_list = data['mjd'][ind_list].tolist()\n\t\twarn_index = []\n\t\tn = 0\n\t\twhile n+1<len(mjd_exit_list):\n\t\t\tdelta = mjd_exit_list[n+1]-mjd_exit_list[n]\n\n\t\t\tif delta>check_delta:\n# \t\t\t\tprint('alone: ',delta) #test\n\t\t\t\twarn_index.append(ind_list[n])\n\t\t\tif n == 0 and delta>260:\n\t\t\t\trecord_row.append(ind_list[0])\n\t\t\tn +=1\n\t\tif len(warn_index)!=0:\n\t\t\tt = 0\n\t\t\twhile t+1<len(warn_index):\n\t\t\t\tif t == 0:\n\t\t\t\t\tif warn_index[0]-ind_list[0]<=min_size: record_row+=list(range(ind_list[0],warn_index[0]+1))\n\n\t\t\t\tif warn_index[t+1]-warn_index[t]<=min_size: \n\t\t\t\t\trecord_row += list(range(warn_index[t]+1, warn_index[t+1]+1))\n\t\t\t\tt +=1\n\t\t\tif len(warn_index)==1:\n\t\t\t\tif warn_index[0]-ind_list[0]<=min_size: record_row+=list(range(ind_list[0],warn_index[0]+1))\n    \n\tif len(record_row)!= 0: print('yes there are alone mjd!')\n\n\tnew_data = data.drop(index=list(set(record_row)), axis=0).reset_index(drop=True)\n\n\treturn new_data\n\n\n\n\ndef combine_narrow_mjd(id_list = list, data=list, check_delta = 0.65):\n\t'''\n\tThis function is used for combine the data points whose mjds are closed. \n\tFor example, 2 data points in one mjd.\n\tInput:\n\t- id_list: the id list of all objects\n\t- data: the data from the preprocessed file\n\t- check_delta: if the difference between two neighboring data points' mjd is smaller than this value, the later one will be removed.\n\tReturns:\n\t- data: the modified data\n\n\t'''\n\trecord_row = []\n\tfor i in id_list:\n\t\tind_list = data[data.object_id == i].index.tolist()\n\t\tmjd_exit_list = data['mjd'][ind_list].tolist()\n\t\tn = 0\n\t\twhile n+1<len(mjd_exit_list):\n\t\t\tdelta = mjd_exit_list[n+1]-mjd_exit_list[n]\n\t\t\tif delta<check_delta:\n\t\t\t\tprint('combine: ', delta)\n\t\t\t\tprint('id: ', i, ' mjd: ', mjd_exit_list, 'mjd_n+1: ', mjd_exit_list[n+1], 'mjd_n: ', mjd_exit_list[n])\n\t\t\t\trecord_row.append(ind_list[n])\n\t\t\tn +=1\n\tnew_data = data.drop(index=list(set(record_row)), axis=0).reset_index(drop=True)\n\treturn new_data\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def suitable_gap_length(data = list):\n    '''\n    This function is used for light curves with non-identical gap boundaries.\n    '''\n    combine = data.groupby(data['object_id'])\n    id_list = []\n    for _id, group in combine:\n        id_list.append(_id)\n    \n    data = remove_alone_mjd(id_list, data, check_delta = 100, min_size = 3)\n#     data = combine_narrow_mjd(id_list, data, check_delta = 0.65)\n    \n    # find suitable gap length\n    \n    # assume the initial max gap length is 300\n    max_gap = 300\n    \n    for i in id_list:\n        ind_list = data[data.object_id == i].index.tolist()\n        mjd_exit_list = data['mjd'][ind_list].tolist()\n        m = 0\n        obj_delta_time = []\n        while m < len(mjd_exit_list)-1:\n            delta = mjd_exit_list[m+1]-mjd_exit_list[m]\n            obj_delta_time.append(delta)\n            m +=1 \n\n        while len([x for x in obj_delta_time if x >= max_gap]) < 2 and max_gap > 0:\n            max_gap -=1\n        print(max_gap)\n#         if max_gap == 7:\n#             print(i)\n\n#     print(max_gap) \n# find it!! about 90 、88!\n            \n            \n            \n    \n    \n    ","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nimport matplotlib.pyplot as plt \nimport random\n\n\nmeta = pd.read_csv('../input/plasticcunblindeddatasets/plasticc_test_metadata.csv')\nwdf_obj = meta[(meta.ddf_bool<1)]['object_id'].tolist()\ndel meta\ndata = pd.read_csv('../input/plasticc-converted-datasets/converted_test_batch2.csv')\ndata = data[(data.object_id.isin(wdf_obj))]\n\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"combine = data.groupby(data['object_id'])\nid_list = []\nfor _id, group in combine:\n    id_list.append(_id)\nsub_num = int(len(id_list)/10)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sub_num","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sub_data = data[(data.object_id.isin(id_list[:sub_num]))]\nsuitable_gap_length(sub_data)\n","execution_count":null,"outputs":[]}],"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat":4,"nbformat_minor":4}