Hệ số tương quan Pearson

* Phân tích tương quan Pearson trong SPSS
* Tác giả: Thực hành kinh tế lượng
* Mô tả: Tính hệ số tương quan Pearson giữa hai biến liên tục

* Mở file dữ liệu
GET FILE='C:\data\correlation_data.sav'.

* Thống kê mô tả ban đầu
DESCRIPTIVES VARIABLES=thu_nhap chi_tieu
  /STATISTICS=MEAN STDDEV MIN MAX.

* Kiểm tra phân phối chuẩn
EXAMINE VARIABLES=thu_nhap chi_tieu
  /PLOT=HISTOGRAM NPPLOT
  /STATISTICS=DESCRIPTIVES
  /CINTERVAL 95.

* Tính hệ số tương quan Pearson
CORRELATIONS
  /VARIABLES=thu_nhap chi_tieu
  /PRINT=TWOTAIL NOSIG
  /MISSING=PAIRWISE.

* Biểu đồ phân tán
GRAPH
  /SCATTERPLOT(BIVAR)=thu_nhap WITH chi_tieu
  /MISSING=LISTWISE.

* Kiểm định ý nghĩa với bootstrap
CORRELATIONS
  /VARIABLES=thu_nhap chi_tieu
  /PRINT=TWOTAIL SIG
  /MISSING=PAIRWISE
  /STATISTICS=DESCRIPTIVES XPROD
  /BOOTSTRAP=YES.

* Phân tích tương quan Pearson trong Stata
* Tác giả: Thực hành kinh tế lượng
* Mô tả: Tính hệ số tương quan Pearson và kiểm định ý nghĩa

* Mở dữ liệu
use "C:\data\correlation_data.dta", clear

* Thống kê mô tả
summarize thu_nhap chi_tieu

* Kiểm tra phân phối chuẩn
sktest thu_nhap chi_tieu
swilk thu_nhap chi_tieu

* Tính hệ số tương quan Pearson
correlate thu_nhap chi_tieu

* Tính hệ số tương quan với kiểm định ý nghĩa
pwcorr thu_nhap chi_tieu, sig

* Ma trận tương quan
corr thu_nhap chi_tieu, means

* Biểu đồ phân tán
scatter chi_tieu thu_nhap || lfit chi_tieu thu_nhap, ///
  title("Mối quan hệ giữa thu nhập và chi tiêu") ///
  xlabel(, grid) ylabel(, grid) ///
  legend(order(1 "Quan sát" 2 "Đường hồi quy"))

* Kiểm định tương quan với bootstrap
bootstrap, reps(1000): correlate thu_nhap chi_tieu

* Xuất kết quả
outreg2 using correlation_results, excel replace

# Phân tích tương quan Pearson trong R
# Tác giả: Thực hành kinh tế lượng
# Mô tả: Tính hệ số tương quan Pearson và visualize kết quả

# Tải các thư viện cần thiết
library(corrplot)      # Vẽ ma trận tương quan
library(ggplot2)       # Vẽ biểu đồ
library(psych)         # Thống kê mô tả
library(nortest)       # Kiểm định phân phối chuẩn
library(Hmisc)         # Tương quan với khoảng tin cậy

# Đọc dữ liệu
data <- read.csv("C:/data/correlation_data.csv")

# Xem cấu trúc dữ liệu
str(data)
head(data)

# Thống kê mô tả
summary(data)
describe(data)

# Kiểm tra phân phối chuẩn
shapiro.test(data$thu_nhap)
shapiro.test(data$chi_tieu)

# Kiểm định Anderson-Darling
ad.test(data$thu_nhap)
ad.test(data$chi_tieu)

# Tính hệ số tương quan Pearson
cor(data$thu_nhap, data$chi_tieu, method = "pearson")

# Kiểm định ý nghĩa
cor.test(data$thu_nhap, data$chi_tieu, method = "pearson")

# Tương quan với khoảng tin cậy
rcorr(as.matrix(data[,c("thu_nhap", "chi_tieu")]), type = "pearson")

# Ma trận tương quan
cor_matrix <- cor(data[,c("thu_nhap", "chi_tieu")])
print(cor_matrix)

# Vẽ ma trận tương quan
corrplot(cor_matrix, method = "circle", type = "upper", 
         order = "hclust", tl.cex = 0.8, tl.col = "black")

# Biểu đồ phân tán với đường hồi quy
ggplot(data, aes(x = thu_nhap, y = chi_tieu)) +
  geom_point(alpha = 0.7, size = 2) +
  geom_smooth(method = "lm", se = TRUE, color = "red") +
  labs(title = "Mối quan hệ giữa thu nhập và chi tiêu",
       x = "Thu nhập (triệu VND)",
       y = "Chi tiêu (triệu VND)") +
  theme_minimal()

# Kiểm tra các giả thiết
# 1. Tính tuyến tính
plot(data$thu_nhap, data$chi_tieu, 
     main = "Biểu đồ phân tán",
     xlab = "Thu nhập", ylab = "Chi tiêu")
abline(lm(chi_tieu ~ thu_nhap, data = data), col = "red")

# 2. Phân phối chuẩn
qqnorm(data$thu_nhap, main = "Q-Q Plot Thu nhập")
qqline(data$thu_nhap)

qqnorm(data$chi_tieu, main = "Q-Q Plot Chi tiêu")
qqline(data$chi_tieu)

# Bootstrap confidence interval
library(boot)
cor_boot <- function(data, indices) {
  sample_data <- data[indices, ]
  return(cor(sample_data$thu_nhap, sample_data$chi_tieu))
}

boot_results <- boot(data, cor_boot, R = 1000)
boot.ci(boot_results, type = "perc")

# Xuất kết quả
write.csv(cor_matrix, "correlation_results.csv")

# Phân tích tương quan Pearson trong Python
# Tác giả: Thực hành kinh tế lượng
# Mô tả: Tính hệ số tương quan Pearson với visualization

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
from scipy import stats
from scipy.stats import pearsonr
import warnings
warnings.filterwarnings('ignore')

# Cài đặt hiển thị tiếng Việt
plt.rcParams['font.family'] = 'DejaVu Sans'
plt.rcParams['figure.figsize'] = (10, 6)

# Đọc dữ liệu
data = pd.read_csv('C:/data/correlation_data.csv')

# Xem cấu trúc dữ liệu
print("Cấu trúc dữ liệu:")
print(data.info())
print("\nThống kê mô tả:")
print(data.describe())

# Kiểm tra missing values
print("\nMissing values:")
print(data.isnull().sum())

# Kiểm tra phân phối chuẩn
def check_normality(data, column):
    """Kiểm tra phân phối chuẩn"""
    statistic, p_value = stats.shapiro(data[column])
    print(f"\nKiểm định Shapiro-Wilk cho {column}:")
    print(f"Statistic: {statistic:.4f}")
    print(f"P-value: {p_value:.4f}")
    if p_value > 0.05:
        print("Dữ liệu có phân phối chuẩn")
    else:
        print("Dữ liệu không có phân phối chuẩn")

check_normality(data, 'thu_nhap')
check_normality(data, 'chi_tieu')

# Tính hệ số tương quan Pearson
correlation, p_value = pearsonr(data['thu_nhap'], data['chi_tieu'])
print(f"\nHệ số tương quan Pearson: {correlation:.4f}")
print(f"P-value: {p_value:.4f}")

# Diễn giải kết quả
def interpret_correlation(r):
    """Diễn giải mức độ tương quan"""
    abs_r = abs(r)
    if abs_r < 0.2:
        return "Rất yếu"
    elif abs_r < 0.4:
        return "Yếu"
    elif abs_r < 0.6:
        return "Đáng kể"
    elif abs_r < 0.8:
        return "Mạnh"
    else:
        return "Rất mạnh"

print(f"Mức độ tương quan: {interpret_correlation(correlation)}")

# Ma trận tương quan
corr_matrix = data[['thu_nhap', 'chi_tieu']].corr()
print("\nMa trận tương quan:")
print(corr_matrix)

# Visualization
# 1. Biểu đồ phân tán
plt.figure(figsize=(12, 8))

plt.subplot(2, 2, 1)
plt.scatter(data['thu_nhap'], data['chi_tieu'], alpha=0.7)
plt.plot(np.unique(data['thu_nhap']), 
         np.poly1d(np.polyfit(data['thu_nhap'], data['chi_tieu'], 1))(np.unique(data['thu_nhap'])), 
         'r--', linewidth=2)
plt.xlabel('Thu nhập (triệu VND)')
plt.ylabel('Chi tiêu (triệu VND)')
plt.title('Biểu đồ phân tán với đường hồi quy')
plt.grid(True, alpha=0.3)

# 2. Heatmap tương quan
plt.subplot(2, 2, 2)
sns.heatmap(corr_matrix, annot=True, cmap='coolwarm', center=0,
            square=True, linewidths=0.5)
plt.title('Ma trận tương quan')

# 3. Histogram
plt.subplot(2, 2, 3)
plt.hist(data['thu_nhap'], bins=20, alpha=0.7, color='skyblue', edgecolor='black')
plt.xlabel('Thu nhập (triệu VND)')
plt.ylabel('Tần số')
plt.title('Phân phối thu nhập')

plt.subplot(2, 2, 4)
plt.hist(data['chi_tieu'], bins=20, alpha=0.7, color='lightgreen', edgecolor='black')
plt.xlabel('Chi tiêu (triệu VND)')
plt.ylabel('Tần số')
plt.title('Phân phối chi tiêu')

plt.tight_layout()
plt.show()

# Q-Q Plot để kiểm tra phân phối chuẩn
fig, axes = plt.subplots(1, 2, figsize=(12, 5))

stats.probplot(data['thu_nhap'], dist="norm", plot=axes[0])
axes[0].set_title('Q-Q Plot Thu nhập')
axes[0].grid(True, alpha=0.3)

stats.probplot(data['chi_tieu'], dist="norm", plot=axes[1])
axes[1].set_title('Q-Q Plot Chi tiêu')
axes[1].grid(True, alpha=0.3)

plt.tight_layout()
plt.show()

# Confidence interval bằng bootstrap
def bootstrap_correlation(x, y, n_bootstrap=1000):
    """Tính khoảng tin cậy bằng bootstrap"""
    bootstrap_corrs = []
    n = len(x)
    
    for _ in range(n_bootstrap):
        indices = np.random.choice(n, size=n, replace=True)
        x_boot = x[indices]
        y_boot = y[indices]
        corr, _ = pearsonr(x_boot, y_boot)
        bootstrap_corrs.append(corr)
    
    return np.array(bootstrap_corrs)

boot_corrs = bootstrap_correlation(data['thu_nhap'].values, data['chi_tieu'].values)
ci_lower = np.percentile(boot_corrs, 2.5)
ci_upper = np.percentile(boot_corrs, 97.5)

print(f"\nKhoảng tin cậy 95% cho hệ số tương quan:")
print(f"[{ci_lower:.4f}, {ci_upper:.4f}]")

# Regression analysis bổ sung
from sklearn.linear_model import LinearRegression
from sklearn.metrics import r2_score

X = data[['thu_nhap']]
y = data['chi_tieu']

model = LinearRegression()
model.fit(X, y)

y_pred = model.predict(X)
r2 = r2_score(y, y_pred)

print(f"\nPhân tích hồi quy bổ sung:")
print(f"R-squared: {r2:.4f}")
print(f"Hệ số góc: {model.coef_[0]:.4f}")
print(f"Hệ số chặn: {model.intercept_:.4f}")

# Xuất kết quả
results = pd.DataFrame({
    'Metric': ['Correlation', 'P-value', 'R-squared', 'CI_Lower', 'CI_Upper'],
    'Value': [correlation, p_value, r2, ci_lower, ci_upper]
})

results.to_csv('correlation_results.csv', index=False)
print("\nKết quả đã được xuất ra file correlation_results.csv")