Разница между PCA в R и Python - PullRequest
Новая
       24

Разница между PCA в R и Python

1 голос
/ 19 февраля 2020

Мой коллега и я сравниваем PCA, выполненный в R и Python. Она использует princomp() функцию в R, тогда как я использую PCA от sklearn.decomposition. Что интересно, после подгонки масштабированных данных в pca.fit и применения уменьшения размеров мы получаем совершенно разные нагрузки для ПК1. Разница в знаке - все компоненты с положительным (отрицательным) знаком в R имеют отрицательный (положительный) знак в Python. Затем мы вписываем результат в модель линейной регрессии и, что неудивительно, мы получаем другой знак коэффициента (тогда как точка пересечения одинакова). У нас также есть то же значение R в квадрате для модели.

Python: 

Coefficients: 
 [-0.17878575]
Intercept: 
 0.6805798284875717

R: 

Coefficients: 
 [0.17878575]
Intercept: 
 0.6805798284875717

Python код выглядит следующим образом: 

scaled_data = pd.DataFrame(preprocessing.scale(df_target_features),columns = df_target_features.columns) 

# mean vector to compute covariance matrix
mean_vector = np.mean(scaled_data, axis=0)

# covarariance matrix, representation of variance of all features in the original data set
cov_matrix = (scaled_data - mean_vector).T.dot((scaled_data - mean_vector)) / (scaled_data.shape[0]-1)

# performing eigen decomposition on covariance matrix to find eigen values and eigen vector to find where is the maximum variance in the data
eigen_vals, eigen_vecs = np.linalg.eig(cov_matrix)

# visually confirm that the list is correctly sorted by decreasing eigenvalues
eigen_pairs = [(np.abs(eigen_vals[i]), eigen_vecs[:, i]) for i in range(len(eigen_vals))]
for i in eigen_pairs:
    print (i[0])




# initialize pca object
components = 1
pca = PCA(n_components=components)


# Fit the model with X and apply the dimensionality reduction on X.The new components are just the two main dimensions of variation.
pca.fit(scaled_data)

# X features - principal components
X_reg = pca.fit_transform(scaled_data)
print(pca.explained_variance_ratio_)

explained_variance = pca.explained_variance_ratio_
pca_components = pca.components_

Y = Y.reset_index(drop=True)

# constructing data frame with principal components (result of dimension reduction)
principal_df = pd.DataFrame(data=X_reg, columns=['PC_' + str(i) for i in range(1, components + 1)])

# we can concatenate principal components with target (Y) variable
final_df = pd.concat([principal_df, Y], axis=1)

comp_df = pd.DataFrame(pca.components_, columns=scaled_data.columns, index=['PC_' + str(i) for i in range(1, components + 1)])

# explained variance

plt.plot(np.cumsum(pca.explained_variance_ratio_))
plt.xlabel('Number of components')
plt.ylabel('Cumulative explained variance')
plt.show()


# Create linear regression object

regr = linear_model.LinearRegression()

# Fit 
regr.fit(X_reg, Y)

# calibration 
y_c = regr.predict(X_reg)

sc = regr.score(X_reg, Y)

# The coefficients
# coefficients = pd.concat([pd.DataFrame(X_reg.columns),pd.DataFrame(np.transpose(regr.coef_))], axis = 1)
print('Coefficients: \n', regr.coef_)
print('Intercept: \n', regr.intercept_)


# Cross-validation; The parameter cv=10 means that we divide the data in 10 parts and hold out one part for cross validation.
y_cv = cross_val_predict(regr, X_reg, Y, cv=10)
# Calculate scores for calibration and cross-validation
score_c = r2_score(Y, y_c)
score_cv = r2_score(Y, y_cv)
# Calculate mean square error for calibration and cross validation
mse_c = mean_squared_error(Y, y_c)
mse_cv = mean_squared_error(Y, y_cv)

Как мы можем решить, у кого есть правильные результаты и когда мне следует сначала провести расследование?

...