Samfunnsøkonomi med datavitenskap

#!/usr/bin/python
# -*- coding: UTF-8 -*-

#Uncomment next line to run in goolge colab
#!pip install geopandas

import numpy as np
import geopandas as gpd
import pandas as pd
import matplotlib.pyplot as plt 
import copy
import urllib
import os
np. set_printoptions(suppress=True)

fp_geodata="https://titlon.uit.no/hht/data/continents.zip"

#This module creates a set of image files that can be used together with 
#a slider tool in a web page to create a graphical representation of population 
#developement the last 2000 years. See the source code of the web page to see
#how to use it in a web page
def main():
    #creating the plot
    fig,ax=plt.subplots()
    
    #creating the world:
    world=world_map()
    
    #defining the years that will be plotted:
    rng=[0, 200, 400, 600, 800, 1000, 1100, 1200, 1300, 1400, 1450, 1500, 1550, 
            1600, 1650, 1700, 1750, 1800, 1825, 1850, 1875, 1900, 1910, 1920, 1930, 1940, 1950, 1960, 
            1970, 1980, 1990, 2000, 2010, 2020, 2030]
    
    #making directory for image files:
    try:
        os.mkdir("./img")
    except:
        pass#just continue if it all ready exist
    
    #iterating over years:
    k=0
    for year in rng:
        w=world.world_calc(year)				#calculating the world for a given year. w is a geopandas map object
        plt.cla()								#clears the plot for a new one
        plt.axis('off')							#do not plot the axes, they make it look less like a map
        plt.axis([-150, 200, -110, 110])		#fix the axis so that the size is constant
        
        #plot the map and save it as a jpg image. This will
        #throw an exception each time, but the figure will be saved fine:
        w.plot('CONTINENT',ax=ax,aspect=0.9)
        fig.savefig(f'./img/fig{k}.jpg',dpi=200,bbox_inches='tight')
        k+=1
    create_html_file(rng)
        
def create_html_file(rng):
    """This function creates a html-file with a slider that select the image for the desired year"""
    n=len(rng)
    s=f"<input type='range' min='0' max='{n-1}' value='34' class='slider' id='myRange'><br><br><br><br>\r\n"
    #iterating to create a html-tag for each image:
    for k in range(len(rng)):#paste the result of this into the html code
        s+=f"<img class='slides' src='./img/fig{k}.jpg'>\r\n"
        
    #building the html-code from a string and two files:
    html=s+f"\r\n</body>\r\n<script>\r\n	var array = JSON.parse('{rng}');"	
    html_fileI = urllib. request.urlopen(r"https://titlon.uit.no/hht/worldI.html").read().decode('utf-8')
    html_fileII = urllib. request.urlopen(r"https://titlon.uit.no/hht/worldII.html").read().decode('utf-8')
    
    #creating the target file and writing to it, before closing:
    result_file=open("world population.html",'w+',encoding='utf-8',newline='')
    result_file.write(html_fileI+html+html_fileII)
    result_file.close()


class world_map:
    """This is a class that stores general properties of the world, 
    and recalcualtes it for a given year"""
    def __init__(self):
        base_year=2030
        
        #read the world geodata:
        world_map_orig=gpd.read_file(fp_geodata)
        world_map_orig=world_map_orig.drop(7)#dropping antarctica
        
        #we need an original and a copy we can play with:
        self.world_map_orig=world_map_orig.drop(5)#dropping oceania
        self.world_map=copy.deepcopy(self.world_map_orig)
        
        n,k=world_pop.shape															#our data has n continents and k years
        self.world_pop=world_pop/get_pop(world_pop, base_year)						#obtaining the population relative to the base_year
        self.world_pop=n*self.world_pop/np.sum(self.world_pop,0).reshape((1,k))		#obtaining the population relative to average
        
    def world_calc(self,year):
        #interpolating:
        p=get_pop(self.world_pop, year)
        
        #iterate though the continents to calcualte the new geometry:
        for i in range(len(names)):
            c=self.world_map_orig[self.world_map['CONTINENT']==names[i]]
            self.world_map['geometry'][self.world_map['CONTINENT']==names[i]]=c.scale(p[i],p[i],p[i],tuple(center[i]))
            
        return self.world_map

def get_pop(pop,year):
    """Calculates the i interpolation of the populations between datapoints"""
    n,k=pop.shape 			#our data has n continents and k years
    
    #ensuring the year lies in the allowed range:
    year=min((year,2030))	
    year=max((year,1))
    
    #if the year is in the data, we don't need interpolation:
    if year in year_pop:
        i=list(year_pop).index(year)
        return pop[:,i:i+1]
    
    #interpolating:
    i=np.nonzero((year>year_pop[:-1])* (year<year_pop[1:]))[0][0]
    a=np.log(pop[:,i:i+1])
    b=np.log(pop[:,i+1:i+2])
    ta=year_pop[i]
    tb=year_pop[i+1]
    pop_est=np.exp(a+((b-a)/(tb-ta))*(year-ta))
    return pop_est


#The data:
year_pop=np.array(
        [1   ,   1000,   1500,   1600,   1700,   1820,   1913,   2000,   2030],dtype=int)
world_pop=np.array([
        [  17. ,   32. ,   47. ,   55. ,   61. ,   74. ,  125. ,  798. , 1449. ],
        [ 336. ,  366. ,  578. ,  758. ,  804. , 1420. , 1956. , 7210. , 9580. ],
        [   0.4,    0.4,    0.6,    0.6,    0.6,    0.4,    5. ,   23. ,   28. ],
        [  34. ,   40. ,   78. ,  112. ,  127. ,  224. ,  498. ,  742. ,  829. ],
        [   1. ,    1. ,    2. ,    2. ,    1. ,   11. ,  105. ,  314. ,  413. ],
        [  27. ,   31. ,   50. ,   20. ,    7. ,   11. ,   38. ,  348. ,  461. ]])

names=['Africa', 'Asia', 'Australia', 'Europe', 'North America', 'South America']


#The center of each continent:
center=np.array([	
            [18.276615, 	6.426117],
            [95.750373,		44.947893],
            [134.489563,	-25.734969],
            [28.24904,		55.854069],
            [-92.320436,	56.515209],
            [-60.792113,	-15.173251]])


main()