import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import statsmodels.api as sm
import seaborn as sns
sns.set() 
##use sklearn
from sklearn.linear_model import LinearRegression   

data = pd.read_csv('C:\\Users\\Python_practice\\1.02. Multiple linear regression.csv')   

data.describe()   ##The dataset has 84 samples

x = data[['SAT','Rand 1,2,3']]
y = data['GPA']

reg = LinearRegression()
reg.fit(x,y)

reg.coef_

reg.intercept_   # We get y = 0.296 + 0.0017*SAT - 0.0083*(Rand 1,2,3)

reg.score(x,y)    ##This is R-square, not adjusted R-square. Usually we use adjusted one to analyze Multilinear 

##If We use feature with little explinatory power, the R-square would increase. 
##Thus we need to penalize this excessive usage through the adjusted R-square 

x.shape  ##n=84(the number of observations), p=2(the number of predictors)

r2 = reg.score(x,y)
n = x.shape[0]
p = x.shape[1]
Rsquare_adj = [1 - (1 - r2)*(n-1)/(n-p-1)]
Rsquare_adj

##Conclusion：Adjusted R-square is considerably less than R-square 
##Thus one or more of the predictors have little or no explinatory power
##We need to eliminate those unnecessary features  
##p-value>0.05, disregard it, in sklearn called f_regression

from sklearn.feature_selection import f_regression

f_regression(x,y)   ##second array is p-value

p_value = f_regression(x,y)[1]
p_value

p_value.round(3)   ##We don't need so many digits, here we take only 3 digits

##Make a conclusion table
reg_summary = pd.DataFrame(data = x.columns.values, columns = ['Features'])
reg_summary['Coefficients'] = reg.coef_ 
reg_summary['p-value'] = p_value.round(3)
reg_summary