07_OptimizeClustersAllObj.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Clustering using all features"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    " %matplotlib inline\n",
    "import seaborn\n",
    "import matplotlib.pyplot as plt\n",
    "from kmeans import kmeans, randomsample\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import random\n",
    "import os\n",
    "from BorealWeights import BorealWeightedProblem\n",
    "from pyomo.opt import SolverFactory\n",
    "from gradutil import *\n",
    "seedn = 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "revenue, carbon, deadwood, ha = init_boreal()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "X = np.concatenate((carbon.values, ha.values, deadwood.values, revenue.values), axis=1)\n",
    "np.shape(X)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "Y = np.concatenate((carbon.values, ha.values, deadwood.values, revenue.values), axis=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.nanmin(ha.values)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.nanmax(ha.values)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.nanmin(X[:,7:14])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%time\n",
    "random.seed(seedn)\n",
    "np.random.seed(seedn)\n",
    "# preprocessing to add penalty for Nan values\n",
    "X[np.isnan(X)] = np.nanmin(X) - np.nanmax(X)\n",
    "randomcenters = randomsample(X, 50)\n",
    "centers, xtoc, dist = kmeans(X,\n",
    "                             randomcenters,\n",
    "                             delta=.00001,\n",
    "                             maxiter=100,\n",
    "                             metric='cosine',\n",
    "                             verbose=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%time\n",
    "C = centers.copy()\n",
    "nvar = len(X)\n",
    "weights = np.array([sum(xtoc==i) for i in range(len(C))])/nvar"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "Ccarbon = C[:,0:7]\n",
    "Cha = C[:,7:14]\n",
    "Cdeadwood = C[:,14:21]\n",
    "Crevenue = C[:,21:]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%time\n",
    "opt = SolverFactory('glpk')\n",
    "\n",
    "clustProblemCarbon = BorealWeightedProblem(Ccarbon,weights,nvar)\n",
    "resCarbon = opt.solve(clustProblemCarbon.model, False)\n",
    "\n",
    "clustProblemHa = BorealWeightedProblem(Cha,weights,nvar)\n",
    "resHA = opt.solve(clustProblemHa.model, False)\n",
    "\n",
    "clustProblemDeadwood = BorealWeightedProblem(Cdeadwood,weights,nvar)\n",
    "resDeadwood = opt.solve(clustProblemDeadwood.model, False)\n",
    "\n",
    "clustProblemRevenue = BorealWeightedProblem(Crevenue,weights,nvar)\n",
    "resRevenue = opt.solve(clustProblemRevenue.model, False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "carbonSurrogateList = res_to_list(clustProblemCarbon.model)\n",
    "haSurrogateList = res_to_list(clustProblemHa.model)\n",
    "deadwoodSurrogateList = res_to_list(clustProblemDeadwood.model)\n",
    "revenueSurrogateList = res_to_list(clustProblemRevenue.model)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "resultSurrogateCarbon = cluster_to_value(Ccarbon, carbonSurrogateList, weights)\n",
    "resultSurrogateHa = cluster_to_value(Cha, haSurrogateList, weights)\n",
    "resultSurrogateDeadwood = cluster_to_value(Cdeadwood, deadwoodSurrogateList, weights)\n",
    "resultSurrogateRev = cluster_to_value(Crevenue, revenueSurrogateList, weights)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Results straight from the surrogate values:')\n",
    "print(\"(i) Harvest revenues {:.0f} M\u20ac\".format(resultSurrogateRev/1000000))\n",
    "print(\"(ii) Carbon storage {:.0f} x 10\u00b3 MgC\".format(resultSurrogateCarbon/1e+3))\n",
    "print(\"(iii) Deadwood index {:.0f} m\u00b3\".format(resultSurrogateDeadwood))\n",
    "print(\"(iv) Combined Habitat {:.0f}\".format(resultSurrogateHa))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Results straight from the surrogate values:')\n",
    "print(\"(i) Harvest revenues {:.0f} M\u20ac\".format(res_value(resRevenue)/1000000))\n",
    "print(\"(ii) Carbon storage {:.0f} x 10\u00b3 MgC\".format(res_value(resCarbon)/1e+3))\n",
    "print(\"(iii) Deadwood index {:.0f} m\u00b3\".format(res_value(resDeadwood)))\n",
    "print(\"(iv) Combined Habitat {:.0f}\".format(res_value(resHA)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "resultOriginCarbon = clusters_to_origin(X[:,:7], xtoc, carbonSurrogateList)\n",
    "resultOriginHa = clusters_to_origin(X[:,7:14], xtoc, haSurrogateList)\n",
    "resultOriginDeadwood = clusters_to_origin(X[:,14:21], xtoc, deadwoodSurrogateList)\n",
    "resultOriginRev = clusters_to_origin(X[:,21:], xtoc, revenueSurrogateList)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Results when surrogate mapped to real values:')\n",
    "print(\"(i) Harvest revenues {:.0f} M\u20ac\".format(resultOriginRev/1000000))\n",
    "print(\"(ii) Carbon storage {:.0f} x 100 MgC\".format(resultOriginCarbon/100))\n",
    "print(\"(iii) Deadwood index {:.0f} m3\".format(resultOriginDeadwood))\n",
    "\n",
    "print(\"(iv) Combined Habitat {:.0f}\".format(resultOriginHa))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "These are the correct results from the original single objective optimizations:\n",
    "- (i) Harvest revenues 250 M\u20ac\n",
    "- (ii) Carbon storage 44490 x 100 MgC\n",
    "- (iii) Deadwood index 218153 m3\n",
    "- (iv) Combined Habitat 20225\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "There is not really so great differences than I thought before, so no worries actually!"
   ]
  }
 ],
 "metadata": {
  "git": {
   "suppress_outputs": true
  },
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.5.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}