<jupyter_start><jupyter_text>Partial Least Square Function<jupyter_code>def Normalize(X):
    '''func to centerize and normalize the dataset,dataset should be numpy array'''
    
    return (X - np.mean(X, axis = 0))/(np.std(X, axis = 0))
def norm(x):
    sum=0
    for i in x:
        sum = sum + i**2
    
    return np.sqrt(sum)
def PLS(X,Y,ncomponents,tol=1e-6):
    E,F=X,Y
    T = []
    W = []
    Q = []
    U = []
    P = []
    B = []
    rY, cY = Y.shape
    rX, cX = X.shape
    for i in range(ncomponents):
        index=np.random.choice(range(Y.shape[1]))
        #u=Y[:,index]
        u=np.random.rand(rY)
        counter = 0
        while(True):
            w = E.T@u
            w = w/norm(w)
            t = E@w
            t = t/norm(t)
            q = F.T@t
            q = q/norm(q)
            u = F@q
            
            if counter==0:
                tlast=t
            elif norm(tlast-t)<tol:
                break
            else:
                tlast=t
            
            counter=counter+1
        
        b = t.T@u
        p = E.T@t
        
        B.append(b)
        T.append(t)
        P.append(p)
        W.append(w)
        Q.append(q)
        U.append(u)
        E = E-t.reshape(-1,1)@p.reshape(1,-1)
        F = F-b*t.reshape(-1,1)@q.reshape(1,-1)
    
        
    return (np.array(T),np.array(P),np.array(W),np.array(Q),np.array(U),np.diag(B))<jupyter_output><empty_output><jupyter_text>Test Function on Wine Data<jupyter_code>#Example1 Data : Wine
X1 = np.array([[7, 7, 13, 7], 
              [4, 3, 14, 7],
              [10, 5, 12, 5],
              [16, 7, 11, 3],
             [13, 3, 10, 3]])

Y1 = np.array([[14, 7, 8], 
              [10, 7, 6],
              [8, 5, 5],
              [2, 4, 7],
              [6, 2, 4]])
X1=Normalize(X1)
Y1=Normalize(Y1)
[T, P, W, Q, U, B] = PLS(X1,Y1,2)
P = P.T
Q = Q.T
P
BPLS = la.pinv(P.T)@B@Q.T
BPLS<jupyter_output><empty_output><jupyter_text>Compare OLS and PLS when there is only one solution<jupyter_code># OLS vs PLS
X_sim = np.random.randn(5, 5)
X_sim
Y_sim = np.random.randn(5,1)
Y_sim
X_sim = Normalize(X_sim)
Y_sim = Normalize(Y_sim)
from sklearn.linear_model import LinearRegression

OLS = LinearRegression()
B_O = OLS.fit(X_sim, Y_sim).coef_.T
B_O
[T, P, W, Q, U, B] = PLS(X_sim,Y_sim,5)

P = P.T
Q = Q.T
B_P = la.pinv(P.T)@B@Q.T
B_P
np.allclose(B_O, B_P)
pls = PLSRegression(n_components=5)

pls.fit(X_sim, Y_sim).coef_<jupyter_output>/Applications/anaconda3/lib/python3.7/site-packages/sklearn/cross_decomposition/pls_.py:292: UserWarning: Y residual constant at iteration 4
  warnings.warn('Y residual constant at iteration %s' % k)<jupyter_text>PLS Application & Comparison<jupyter_code>#Import cars data
df = pd.read_excel("/Users/rachelxing/Desktop/STA663/cars_pls_regression.xls")
df.head()
X = df.iloc[:,:-3].to_numpy()
Y = df.iloc[:, -3:].to_numpy()
X.shape, Y.shape
#normalize X and Y
X = Normalize(X)
Y = Normalize(Y)
#PLS + leave one out (20 fold)
kf = KFold(n_splits=20, random_state=None, shuffle=False)
kf.get_n_splits(X)

y_predict_pls = []
y_test_pls = []

for train_index, test_index in kf.split(X):
    
    X_train, X_test = X[train_index], X[test_index]
    y_train, y_test = Y[train_index], Y[test_index]
    
    [T, P, W, Q, U, B] = PLS(X_train,y_train,7)
    
    P = P.T
    Q = Q.T
    BPLS = la.pinv(P.T)@B@Q.T
    
    y_test_pls.append(y_test)
    y_predict_pls.append(X_test@BPLS)
y_predict_pls = np.array(y_predict_pls).reshape(20,3)
y_test_pls = np.array(y_test_pls).reshape(20,3)
mean_squared_error(y_test_pls, y_predict_pls)
#OLS on cars data + leave one out 
y_predict_ols = []
y_test_ols = []


for train_index, test_index in kf.split(X):
    
    X_train, X_test = X[train_index], X[test_index]
    y_train, y_test = Y[train_index], Y[test_index]
    
    reg1 = LinearRegression().fit(X_train, y_train[:,0])
    reg2 = LinearRegression().fit(X_train, y_train[:,1])
    reg3 = LinearRegression().fit(X_train, y_train[:,2])
    
    p1 = reg1.predict(X_test)
    p2 = reg2.predict(X_test)
    p3 = reg3.predict(X_test)
    
    y_test_ols.append(y_test)
    y_predict_ols.append([p1 ,p2, p3])



y_predict_ols = np.array(y_predict_ols).reshape(20,3)
y_test_ols = np.array(y_test_ols).reshape(20,3)
mean_squared_error(y_test_ols, y_predict_ols)
#Ridge Regression

#Select best parameter alpha
ridge = Ridge()
parameters = {'alpha' : [1e-10, 1e-8, 1e-4, 1e-3, 1e-2, 1, 5, 10, 20]}
ridge_reg = GridSearchCV(ridge, parameters, scoring = 'neg_mean_squared_error', cv = 20)
ridge_reg.fit(X, Y)
print(ridge_reg.best_params_)
print(ridge_reg.best_score_)

#Ridge Regression

y_predict_ridge = []
y_test_ridge = []

for train_index, test_index in kf.split(X):
    
    X_train, X_test = X[train_index], X[test_index]
    y_train, y_test = Y[train_index], Y[test_index]
    
    reg = Ridge(alpha=5)
    reg.fit(X_train, y_train)
    
    
   
    
    y_test_ridge.append(y_test)
    y_predict_ridge.append(reg.predict(X_test))


y_predict_ridge = np.array(y_predict_ridge).reshape(20,3)
y_test_ridge = np.array(y_test_ridge).reshape(20,3)
mean_squared_error(y_test_ridge, y_predict_ridge)
#Principal Component Regression 
pca = PCA(n_components=7)
pca.fit(X.T)
print(pca.explained_variance_ratio_)
Z = pca.components_.T
X.shape, pca.components_.T.shape
#Regress on Principal components

y_predict_pcr = []
y_test_pcr = []


for train_index, test_index in kf.split(Z):
    
    X_train, X_test = Z[train_index], Z[test_index]
    y_train, y_test = Y[train_index], Y[test_index]
    
    reg1 = LinearRegression().fit(X_train, y_train[:,0])
    reg2 = LinearRegression().fit(X_train, y_train[:,1])
    reg3 = LinearRegression().fit(X_train, y_train[:,2])
    
    p1 = reg1.predict(X_test)
    p2 = reg2.predict(X_test)
    p3 = reg3.predict(X_test)
    
    y_test_pcr.append(y_test)
    y_predict_pcr.append([p1 ,p2, p3])

    
y_predict_pcr = np.array(y_predict_pcr).reshape(20,3)
y_test_pcr = np.array(y_test_pcr).reshape(20,3)

mean_squared_error(y_test_pcr, y_predict_pcr)<jupyter_output><empty_output><jupyter_text>Visualization<jupyter_code>df_test =pd.DataFrame(Y, columns=['N_conscity', 'N_price', 'N_symboling'] )
df_test.head()
df_test[['PLS_conscity', 'PLS_price', 'PLS_symboling']] = pd.DataFrame(y_predict_pls)
df_test.head()
fig, axs = plt.subplots(1,3, figsize = (15, 5))
fig.suptitle('PLS Performance', fontsize=20)
axs[0].scatter(df_test["N_conscity"], df_test["PLS_conscity"] , c = 'black')
axs[0].plot([0, 1], [0, 1], transform=axs[0].transAxes, c = 'black', linestyle='dashed')
axs[0].set_xlabel('Conscity (test)')
axs[0].set_ylabel('Conscity (predict)')



axs[1].scatter(df_test["N_price"], df_test["PLS_price"] , c = 'black')
axs[1].plot([0, 1], [0, 1], transform=axs[1].transAxes, c = 'black', linestyle='dashed')
axs[1].set_xlabel('Price (test)')
axs[1].set_ylabel('Price (predict)')


axs[2].scatter(df_test["N_symboling"], df_test["PLS_symboling"] , c = 'black')
axs[2].plot([0, 1], [0, 1], transform=axs[2].transAxes, c = 'black', linestyle='dashed')
axs[2].set_xlabel('Symboling (test)')
axs[2].set_ylabel('Symboling (predict)')
df_test[['OLS_conscity', 'OLS_price', 'OLS_symboling']] = pd.DataFrame(y_predict_ols)
df_test.head()
fig, axs = plt.subplots(1,3, figsize = (15, 5))
fig.suptitle('OLS Performance', fontsize=20)
axs[0].scatter(df_test["N_conscity"], df_test["OLS_conscity"])
axs[0].plot([0, 1], [0, 1], transform=axs[0].transAxes, linestyle='dashed')
axs[0].set_xlabel('Conscity (test)')
axs[0].set_ylabel('Conscity (predict)')


axs[1].scatter(df_test["N_price"], df_test["OLS_price"] )
axs[1].plot([0, 1], [0, 1], transform=axs[1].transAxes,  linestyle='dashed')
axs[1].set_xlabel('Price (test)')
axs[1].set_ylabel('Price (predict)')


axs[2].scatter(df_test["N_symboling"], df_test["OLS_symboling"] )
axs[2].plot([0, 1], [0, 1], transform=axs[2].transAxes,  linestyle='dashed')
axs[2].set_xlabel('Symboling (test)')
axs[2].set_ylabel('Symboling (predict)')
df_test[['Ridge_conscity', 'Ridge_price', 'Ridge_symboling']] = pd.DataFrame(y_predict_ridge)
df_test.head()
fig, axs = plt.subplots(1,3, figsize = (15, 5))
fig.suptitle('Ridge Performance', fontsize=20)
axs[0].scatter(df_test["N_conscity"], df_test["Ridge_conscity"], c = 'orange' )
axs[0].plot([0, 1], [0, 1], transform=axs[0].transAxes, c = 'orange', linestyle='dashed')
axs[0].set_xlabel('Conscity (test)')
axs[0].set_ylabel('Conscity (predict)')


axs[1].scatter(df_test["N_price"], df_test["Ridge_price"], c = 'orange' )
axs[1].plot([0, 1], [0, 1], transform=axs[1].transAxes, c = 'orange', linestyle='dashed')
axs[1].set_xlabel('Price (test)')
axs[1].set_ylabel('Price (predict)')


axs[2].scatter(df_test["N_symboling"], df_test["Ridge_symboling"], c = 'orange' )
axs[2].plot([0, 1], [0, 1], transform=axs[2].transAxes, c = 'orange', linestyle='dashed')
axs[2].set_xlabel('Symboling (test)')
axs[2].set_ylabel('Symboling (predict)')
df_test[['PCR_conscity', 'PCR_price', 'PCR_symboling']] = pd.DataFrame(y_predict_pcr)
df_test.head()
fig, axs = plt.subplots(1,3, figsize = (15, 5))
fig.suptitle('PCR Performance', fontsize=20)
axs[0].scatter(df_test["N_conscity"], df_test["PCR_conscity"], c = 'navy' )
axs[0].plot([0, 1], [0, 1], transform=axs[0].transAxes, c = 'navy', linestyle='dashed')
axs[0].set_xlabel('Conscity (test)')
axs[0].set_ylabel('Conscity (predict)')


axs[1].scatter(df_test["N_price"], df_test["PCR_price"], c = 'navy' )
axs[1].plot([0, 1], [0, 1], transform=axs[1].transAxes, c = 'navy', linestyle='dashed')
axs[1].set_xlabel('Price (test)')
axs[1].set_ylabel('Price (predict)')


axs[2].scatter(df_test["N_symboling"], df_test["PCR_symboling"], c = 'navy' )
axs[2].plot([0, 1], [0, 1], transform=axs[2].transAxes, c = 'navy', linestyle='dashed')
axs[2].set_xlabel('Symboling (test)')
axs[2].set_ylabel('Symboling (predict)')<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.18122581958867848, 'fraction_numerical': 0.03125636326613724, 'mean_word_length': 4.116145833333333, 'pattern_counts': {'":': 0, '<': 21, '<?xml version=': 0, '>': 20, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '5', 'dir': 'jupyter-structured-clean-dedup', 'id': '9739', 'license': 'MIT', 'n_tokens_mistral': 4438, 'n_tokens_neox': 4244, 'n_words': 697, 'path': 'PLS_Final .ipynb', 'repo_name': 'Juncheng-Dong/Partial-Least-Square', 'size': '9822'}