{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.10.13","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":81000,"databundleVersionId":8812083,"sourceType":"competition"}],"dockerImageVersionId":30746,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"Hey fellow FutureCrop enthusiasts 👋\n\nI now also want to share some code.\n\n## Crop Stages\n\nSome real agronomists in my institute told me that crop development is commonly divided into phenological *growth stages* (see e.g. [this](https://ndawn.ndsu.nodak.edu/help-wheat-growing-degree-days.html) webpage by the North Dakota State University). In these stages plants do different biological things and thus have different needs (like they die from different temperatures and stuff).\n\nCommonly, these growth stages are estimated using [growing degree days (GDD)](https://en.wikipedia.org/wiki/Growing_degree-day). When the GDD value reaches a certain threshold, the next growth stage is predicted to start. This seems like super simple 80s shit, but apparently it is really helpful also in farming practice and aligns quite well with visible phenological traits.\n\nNow it seems useful to utilize this domain knowledge for our yield prediction goal. For example by defining the thresholds and calculating summary statistics over these stages. This approach seems to work really well (thanks a lot to [Monique](https://www.kaggle.com/code/picsoflily/r-script01-load-save-features) 💚)!\n\nAnother way would be to just put the GDD into a 1D conv network, so the model can figure out what to do with it itself. I saw this was done in Ron's nice [CNN example](https://www.kaggle.com/code/ronvbree/cnn-example). I had a similar approach and with a lot of tweaking, I could get a leaderboard score of up to `0.816` with it (good, but not as good as Monique 👿 so maybe defined stages are better? or I just suck at some other aspect of regression 🥲).\n\n## Crop Stage Model\nI had problems with overfitting and thus came back to the idea of defined growth stages. Could a model do better than hard-coded thresholds from 80s literature?\n\nTo try that, I put together a simple learnable stage indicator model (in torch). The input is a timeseries of (cumulative) GDD. The output is `num_stages` timeseries of stage indicators (0 to 1). In order to be trainable using gradient-based optimization, the stage boundaries are smooth (sigmoid). Note that the stages do not end, as I found it unnecessarily complex to implement and it can be achieved by linear combination in a later stage anyways. The parameters are just the thresholds for the stages (in the implementation the `bias` params of the conv layer).\n\nI didn't thoroughly test it, but it seems like using the stage indicators as an input to a regular 1D conv net instead of the raw GDD improved the performance. Maybe someone likes to try and tell me, if it improves performance?\n\n(It would be cool to optimize the GDD params (base temperature and maximum temperature) as well, but I couldn't think of a way to do so in a differentiable way 😐️)","metadata":{}},{"cell_type":"code","source":"import torch\nimport matplotlib.pyplot as plt","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2024-08-21T17:39:04.862258Z","iopub.execute_input":"2024-08-21T17:39:04.863646Z","iopub.status.idle":"2024-08-21T17:39:04.869458Z","shell.execute_reply.started":"2024-08-21T17:39:04.863604Z","shell.execute_reply":"2024-08-21T17:39:04.867740Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"class StageIndicatorModel(torch.nn.Module):\n    def __init__(self, num_stages, sharpness=20):\n        super(StageIndicatorModel, self).__init__()\n        self.conv_start = torch.nn.Conv1d(1, num_stages, kernel_size=1)\n        self.conv_start.weight = torch.nn.Parameter(torch.ones(num_stages, dtype=torch.float32)[:, None, None])\n        self.conv_start.weight.requires_grad = False\n        \n        bias = (torch.arange(num_stages, dtype=torch.float32) - num_stages/2) / (num_stages/4)\n        self.conv_start.bias = torch.nn.Parameter(bias)\n        self.conv_start\n        self.activation = torch.nn.Sigmoid()\n        self.sharpness = sharpness\n        \n    def forward(self, gdd):\n        conv_out = self.conv_start(gdd)\n        conv_out *= self.sharpness\n        \n        out = self.activation(conv_out)\n        \n        return out","metadata":{"execution":{"iopub.status.busy":"2024-08-21T17:39:04.871934Z","iopub.execute_input":"2024-08-21T17:39:04.872403Z","iopub.status.idle":"2024-08-21T17:39:04.889033Z","shell.execute_reply.started":"2024-08-21T17:39:04.872370Z","shell.execute_reply":"2024-08-21T17:39:04.887366Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sim = StageIndicatorModel(\n    num_stages=5,\n    sharpness=10,\n)","metadata":{"execution":{"iopub.status.busy":"2024-08-21T17:39:04.890670Z","iopub.execute_input":"2024-08-21T17:39:04.891041Z","iopub.status.idle":"2024-08-21T17:39:04.906570Z","shell.execute_reply.started":"2024-08-21T17:39:04.891012Z","shell.execute_reply":"2024-08-21T17:39:04.904704Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"This is just an example for a standardized (mean=0, std=1) GDD timeseries.","metadata":{}},{"cell_type":"code","source":"gdd = torch.arange(210, dtype=torch.float32).view(1, -1)\ngdd = (gdd-gdd.mean()) / gdd.std()","metadata":{"execution":{"iopub.status.busy":"2024-08-21T17:39:04.908531Z","iopub.execute_input":"2024-08-21T17:39:04.909084Z","iopub.status.idle":"2024-08-21T17:39:04.921966Z","shell.execute_reply.started":"2024-08-21T17:39:04.909033Z","shell.execute_reply":"2024-08-21T17:39:04.919902Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Looks cool IMHO 😋","metadata":{}},{"cell_type":"code","source":"fig, ax = plt.subplots(1, figsize=(12, 5))\nwith torch.no_grad():\n    ax.imshow(sim(gdd), aspect=\"auto\")","metadata":{"execution":{"iopub.status.busy":"2024-08-21T17:39:04.924921Z","iopub.execute_input":"2024-08-21T17:39:04.925564Z","iopub.status.idle":"2024-08-21T17:39:05.284169Z","shell.execute_reply.started":"2024-08-21T17:39:04.925520Z","shell.execute_reply":"2024-08-21T17:39:05.282809Z"},"trusted":true},"execution_count":null,"outputs":[]}]}