09_CompareClusteringsNclust.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Compare different clusterings with different parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "from gradutil import *\n",
    "from BorealWeights import BorealWeightedProblem\n",
    "from pyomo.opt import SolverFactory"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "seedn = 4\n",
    "opt = SolverFactory('glpk')\n",
    "solutions = ideal()\n",
    "revenue, carbon, deadwood, ha = init_boreal()\n",
    "x = np.concatenate((revenue, carbon, deadwood, ha), axis=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Critical point: what features to use in clustering?\n",
    "- cluster_x in this\n",
    "  - assing Nan:s to be (min - max) so that they separate clearly from everything else\n",
    "  - then scale all the variables(columns) from 0 to 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "norm_revenue = normalize(revenue.values)\n",
    "norm_carbon = normalize(carbon.values)\n",
    "norm_deadwood = normalize(deadwood.values)\n",
    "norm_ha = normalize(ha.values)\n",
    "cluster_x = np.concatenate((norm_revenue, norm_carbon, norm_deadwood, norm_ha), axis=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Critical point: how to normalize data for optimization?\n",
    "- no_nan_x in this  \n",
    "  - assign Nan:s to be as (min - max) per column, so that they significanty different to all others"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "no_nan_x = x.copy()\n",
    "inds = np.where(np.isnan(no_nan_x))\n",
    "no_nan_x[inds] = np.take(np.nanmin(no_nan_x, axis=0) - np.nanmax(no_nan_x, axis=0), inds[1])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Critical point: Number of clusters?\n",
    "- nclust in this"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%time\n",
    "nclust = 5\n",
    "optim_revenue5, optim_carbon5, optim_deadwood5, optim_ha5 = cNopt(x, cluster_x, no_nan_x, opt, nclust, seedn)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Relative differences to original values, 5 clusters')\n",
    "print(\"(i) Harvest revenues difference {:.3f}\".format((optim_revenue5-solutions['revenue'])/solutions['revenue']))\n",
    "print(\"(ii) Carbon storage {:.3f}\".format((optim_carbon5-solutions['carbon'])/solutions['carbon']))\n",
    "print(\"(iii) Deadwood index {:.3f}\".format((optim_deadwood5-solutions['deadwood'])/solutions['deadwood']))\n",
    "print(\"(iv) Combined Habitat {:.3f}\".format((optim_ha5-solutions['ha'])/solutions['ha']))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%time\n",
    "nclust = 25\n",
    "optim_revenue25, optim_carbon25, optim_deadwood25, optim_ha25 = cNopt(x, cluster_x, no_nan_x, opt, nclust, seedn)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Relative differences to original values, 25 clusters')\n",
    "print(\"(i) Harvest revenues difference {:.3f}\".format((optim_revenue25-solutions['revenue'])/solutions['revenue']))\n",
    "print(\"(ii) Carbon storage {:.3f}\".format((optim_carbon25-solutions['carbon'])/solutions['carbon']))\n",
    "print(\"(iii) Deadwood index {:.3f}\".format((optim_deadwood25-solutions['deadwood'])/solutions['deadwood']))\n",
    "print(\"(iv) Combined Habitat {:.3f}\".format((optim_ha25-solutions['ha'])/solutions['ha']))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%time\n",
    "nclust = 50\n",
    "optim_revenue50, optim_carbon50, optim_deadwood50, optim_ha50 = cNopt(x, cluster_x, no_nan_x, opt, nclust, seedn)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Relative differences to original values, 50 clusters')\n",
    "print(\"(i) Harvest revenues difference {:.3f}\".format((optim_revenue50-solutions['revenue'])/solutions['revenue']))\n",
    "print(\"(ii) Carbon storage {:.3f}\".format((optim_carbon50-solutions['carbon'])/solutions['carbon']))\n",
    "print(\"(iii) Deadwood index {:.3f}\".format((optim_deadwood50-solutions['deadwood'])/solutions['deadwood']))\n",
    "print(\"(iv) Combined Habitat {:.3f}\".format((optim_ha50-solutions['ha'])/solutions['ha']))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%time\n",
    "nclust = 100\n",
    "opt = SolverFactory('glpk')\n",
    "optim_revenue100, optim_carbon100, optim_deadwood100, optim_ha100 = cNopt(x, cluster_x, no_nan_x, opt, nclust, seedn)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Relative differences to original values, 100 clusters')\n",
    "print(\"(i) Harvest revenues difference {:.3f}\".format((optim_revenue100-solutions['revenue'])/solutions['revenue']))\n",
    "print(\"(ii) Carbon storage {:.3f}\".format((optim_carbon100-solutions['carbon'])/solutions['carbon']))\n",
    "print(\"(iii) Deadwood index {:.3f}\".format((optim_deadwood100-solutions['deadwood'])/solutions['deadwood']))\n",
    "print(\"(iv) Combined Habitat {:.3f}\".format((optim_ha100-solutions['ha'])/solutions['ha']))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%time\n",
    "nclust = 500\n",
    "opt = SolverFactory('glpk')\n",
    "optim_revenue500, optim_carbon500, optim_deadwood500, optim_ha500 = cNopt(x, cluster_x, no_nan_x, opt, nclust, seedn)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Relative differences to original values, 100 clusters')\n",
    "print(\"(i) Harvest revenues difference {:.3f}\".format((optim_revenue500-solutions['revenue'])/solutions['revenue']))\n",
    "print(\"(ii) Carbon storage {:.3f}\".format((optim_carbon500-solutions['carbon'])/solutions['carbon']))\n",
    "print(\"(iii) Deadwood index {:.3f}\".format((optim_deadwood500-solutions['deadwood'])/solutions['deadwood']))\n",
    "print(\"(iv) Combined Habitat {:.3f}\".format((optim_ha500-solutions['ha'])/solutions['ha']))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As we can see, the solutions get more accurate, when the number of clusters increases. Everything is then as it should be."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Same without Nan:s - Not A Problem Anymore"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Because clustering is causing troubles, lets do the same actions using data without Nan-values."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "z = np.concatenate((revenue.dropna(axis=0, how='any'), \n",
    "                    carbon.dropna(axis=0, how='any'), \n",
    "                    deadwood.dropna(axis=0, how='any'), \n",
    "                    ha.dropna(axis=0, how='any')), axis=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "cluster_z = normalize(z)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%time\n",
    "nclust = 50\n",
    "optim_revenuez5, optim_carbonz5, optim_deadwoodz5, optim_haz5 = cNopt(z, cluster_z, z, opt, nclust, seedn)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "sol_revenue, sol_carbon, sol_deadwood, sol_ha = optimize_all(z, np.ones(len(z)), opt)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "revenue_list = values_to_list(sol_revenue[0].model, z[:,:7])\n",
    "carbon_list = values_to_list(sol_carbon[0].model, z[:,7:14])\n",
    "deadwood_list = values_to_list(sol_deadwood[0].model, z[:,14:21])\n",
    "ha_list = values_to_list(sol_ha[0].model, z[:,21:])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "print('Relative differences to original values (calculated without Nan:s), 50 clusters')\n",
    "print(\"(i) Harvest revenues difference {:.3f}\".format((optim_revenuez5-sum(revenue_list))/sum(revenue_list)))\n",
    "print(\"(ii) Carbon storage {:.3f}\".format((optim_carbonz5-sum(carbon_list))/sum(carbon_list)))\n",
    "print(\"(iii) Deadwood index {:.3f}\".format((optim_deadwoodz5-sum(deadwood_list))/sum(deadwood_list)))\n",
    "print(\"(iv) Combined Habitat {:.3f}\".format((optim_haz5-sum(ha_list))/sum(ha_list)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now all the differences are much smaller. So the problem is still in the way the nan-values are normalized. "
   ]
  }
 ],
 "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
}