Daily_RA_Temp_2019.py

# coding: utf-8

# In[3]:


import csv
import matplotlib.pyplot as plt
import numpy as np
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import PolynomialFeatures
from sklearn.model_selection import cross_val_score
from sklearn import datasets, linear_model
from sklearn.metrics import mean_squared_error
import pandas as pd


# In[4]:


filename = "1985_2019.csv"
df = pd.read_csv(filename)

RA=[]

temp = df['RA'][df['MO'] == 6]
RA.append(np.asarray(temp))

temp = df['RA'][df['MO'] == 7]
RA.append(np.asarray(temp))

temp = df['RA'][df['MO'] == 8]
RA.append(np.asarray(temp))

temp = df['RA'][df['MO'] == 9]
RA.append(np.asarray(temp))

#print(RA[0][0])


# In[75]:


#Validation of Model

def init_list_of_objects(size):
    list_of_objects = list()
    for i in range(0,size):
        list_of_objects.append( list() ) #different object reference each time
    return list_of_objects

Y_pred_RA=[None]*4

for I in range(4):
    
#print(ra[1])
    print()
    if I==0:
        print("June: ")
    elif I==1:
        print("July: ")
    elif I==2:
        print("August: ")
    else:
        print("September: ")
        
    day_length=0
    if I==0 or I==3:
        day_length = 30
    else:
        day_length = 31
        
    ra = init_list_of_objects(day_length)
    i=0
    while i<len(RA[I]):
        temp=[]
        temp.append(RA[I][i])
        ra[i%(day_length)].append(temp)
        i=i+1

    #ra[0] for partilcular I indicates day 1 of Ith(starting from June) month ra data from 1985 to 2018
    
    #print(ra[0])
    #print(len(ra[0]))

    pred = []
    for d in range(day_length):

        x_train=[]
        for i in range(32):
            temp = []
            temp.append(i+1)
            x_train.append(temp)
        #print(x_train) 

        x_test=[[33],[34]]

        y_train = ra[d][:-2]
        
        y_test = ra[d][-2:]

        regr = linear_model.LinearRegression()

        regr.fit(x_train, y_train)
        
        pred = regr.predict(x_test)
        
        print('day %s' %str(d+1))
        print('Prediction: %.2f' %pred[0][0] + ' %.2f' %pred[1][0])
        print('True value: %.2f' %y_test[0][0] + ' %.2f' %y_test[1][0])
        print("Mean squared error: %.2f" % mean_squared_error(y_test,pred))
        print()


# In[35]:


#Graph of Validation and Trained Model

def init_list_of_objects(size):
    list_of_objects = list()
    for i in range(0,size):
        list_of_objects.append( list() ) #different object reference each time
    return list_of_objects

Y_pred_RA=[None]*4

for I in range(1):
    
    print()
    
    s=''
    if I==0:
        s="June"
        print("June: ")
    elif I==1:
        s="July"
        print("July: ")
    elif I==2:
        s="August"
        print("August: ")
    else:
        s="September"
        print("September: ")
        
    day_length=0
    if I==0 or I==3:
        day_length = 30
    else:
        day_length = 31
        
    ra = init_list_of_objects(day_length)
    i=0
    while i<len(RA[I]):
        temp=[]
        temp.append(RA[I][i])
        ra[i%(day_length)].append(temp)
        i=i+1

    #ra[0] for partilcular I indicates day 1 of Ith(starting from June) month ra data from 1985 to 2018
    
    #print(ra[0])
    #print(len(ra[0]))

    pred = []
    for d in range(day_length):
        
        g1x_train=[]
        for i in range(31):
            temp = []
            temp.append(i+1985)
            g1x_train.append(temp)

        g1x_test=[[2016],[2017],[2018]]

        g1y_train = ra[d][:-3]

        g1y_test = ra[d][-3:]

        regr = linear_model.LinearRegression()

        regr.fit(g1x_train, g1y_train)

        g1y_pred = regr.predict(g1x_test)

        plt.scatter(g1x_test, g1y_test, color='black')
        plt.plot(g1x_test, g1y_pred, color='red', linewidth=3)
        plt.ylim(top=15,bottom=6)
        plt.xticks([2016,2017,2018])
        plt.legend(('predicted line','original data'))
        plt.title('validation of '+s+' day '+str(d+1))
        plt.show()
        #plt.savefig('graph_after_train.png')
        
        #After validation by testing we use it to predict further year result
        
        g2x_train = []
        for i in range(34):
            temp = []
            temp.append(i+1985)
            g2x_train.append(temp)
        
        g2y_train = ra[d]
        regr.fit(g2x_train,g2y_train)
        g2y_pred = regr.predict(g2x_train)
  
        plt.scatter(g2x_train, g2y_train, color='black')
        plt.plot(g2x_train, g2y_pred, color='blue', linewidth=3)
        plt.ylim(top=15,bottom=6)
        #plt.xticks([2016,2017,2018])
        plt.legend(('predicted line','original data'))
        plt.title('trained graph '+s+' day '+str(d+1))
        plt.show()


# In[61]:


#Prediction for 2019 data daywise

def init_list_of_objects(size):
    list_of_objects = list()
    for i in range(0,size):
        list_of_objects.append( list() ) #different object reference each time
    return list_of_objects

Y_pred_RA=[None]*4

for I in range(4):

    day_length=0
    if I==0 or I==3:
        day_length = 30
    else:
        day_length = 31
        
    ra = init_list_of_objects(day_length)
    i=0
    while i<len(RA[I]):
        temp=[]
        temp.append(RA[I][i])
        ra[i%(day_length)].append(temp)
        i=i+1

    #ra[0] for partilcular I indicates day 1 of Ith(starting from June) month ra data from 1985 to 2018
    
    #print(ra[0])
    #print(len(ra[0]))

    pred = []
    for d in range(day_length):

        x_train=[]
        for i in range(34):
            temp = []
            temp.append(i+1985)
            x_train.append(temp)
        #print(x_train) 

        x_test=[[2019]]

        y_train = ra[d]
        
        y_test = ra[d]

        regr = linear_model.LinearRegression()

        regr.fit(x_train, y_train)
        
        temp_pred = regr.predict(x_test)
    
        pred.append(temp_pred[0][0])
        
    #print(pred)

    Y_pred_RA[I]=pred
    
    print(Y_pred_RA[I])
    print()

    #Y_pred[i][j] = RA of jth day of ith month-------i=0 means june and i=3 means september-----j is till 30 or 31 depends on month

    #print(len(Y_pred_RA[0]))

f = open("RA of 2019 june to sep.txt","w+");

for i in range(4):
    
    day_length
    if i==0 or i==3:
        day_length=30
    else:
        day_length=31
        
    for j in range(day_length):
        f.write("%.2f "%Y_pred_RA[i][j])
        
    f.write("\n")
    
f.close()


# In[4]:


Temp = []

temp = df['TS'][df['MO'] == 6]
Temp.append(np.asarray(temp))

temp = df['TS'][df['MO'] == 7]
Temp.append(np.asarray(temp))

temp = df['TS'][df['MO'] == 8]
Temp.append(np.asarray(temp))

temp = df['TS'][df['MO'] == 9]
Temp.append(np.asarray(temp))


# In[5]:


#Prediction of temperature

Y_pred_Temp=[None]*4

for I in range(4):
    
#print(ra[1])
    day_length=0
    if I==0 or I==3:
        day_length = 30
    else:
        day_length = 31
        
    ts = init_list_of_objects(day_length)
    i=0
    while i<len(Temp[I]):
        temp=[]
        temp.append(Temp[I][i])
        ts[i%(day_length)].append(temp)
        i=i+1

    #print(ra)
    #print(len(ra[0]))

    pred = []
    for d in range(day_length):

        x_train=[]
        for i in range(34):
            temp = []
            temp.append(i+1)
            x_train.append(temp)
        #print(x_train) 

        x_test=[[35]]

        y_train = ts[d];

        regr = linear_model.LinearRegression()

        regr.fit(x_train, y_train)
        
        temp_pred = regr.predict(x_test)
        
        pred.append(temp_pred[0][0])
    #print(pred)

    Y_pred_Temp[I]=pred

#Y_pred[i][j] = RA/temperatue of jth day of ith month-------i=0 means june and i=4 means september-----j is till 30 or 31 depends on month

#print(Y_pred_Temp[0])

f = open("Temp_2019_june_sep.txt","w+");

for i in range(4):
    
    day_length
    if i==0 or i==3:
        day_length=30
    else:
        day_length=31
        
    for j in range(day_length):
        f.write("%f "%Y_pred_Temp[i][j])
        
    f.write("\n")
    
f.close()


# In[25]:


#NDVI
x=[0.548,0.618,0.612,0.551,0.765,0.707,0.726,0.689,0.613,0.669,0.654,0.693,0.451,0.664,0.63,0.128,0.579,0.416,0.682,0.666,0.521,0.67,0.681,0.616]

print(len(x))

avg_x=[]
temp=[]
temp.append(x[0])
for i in range(1,len(x)):
    if i%4==0:
        avg_x.append(sum(temp)/4)
        temp=[]
        temp.append(x[i])
    else:
        temp.append(x[i])

avg_x.append(sum(temp)/4)
print(avg_x)

regr.fit([[1],[2],[3],[4],[5],[6]], avg_x)
        
temp_pred = regr.predict([[7]])

print(temp_pred)
        
#avg_june = [0.2835, 0.24575, 0.21325, 0.2475, 0.27425, 0.172]

#july_ndvi = []