Đồ thị histogram

* Tạo histogram trong SPSS bằng Syntax
* Tác giả: Thực hành kinh tế lượng
* Ví dụ: Phân tích phân phối tuổi lao động

* Nhập dữ liệu mẫu về tuổi
DATA LIST FREE / age.
BEGIN DATA
36 25 38 46 55 68 72 55 36 38
67 45 22 48 91 46 52 61 58 55
END DATA.

* Tạo histogram cơ bản
GRAPH /HISTOGRAM=age.

* Tạo histogram với đường cong phân phối chuẩn
GRAPH /HISTOGRAM(NORMAL)=age.

* Tạo histogram với số bins tùy chỉnh
GRAPH /HISTOGRAM(NORMAL)=age
  /HISTOGRAM FORMAT=INTERVAL WIDTH(AUTO) SCALE(FREQ)
  /TITLE='Phân phối tuổi lao động Việt Nam'.

* Thống kê mô tả kèm theo
FREQUENCIES VARIABLES=age
  /HISTOGRAM NORMAL
  /STATISTICS=STDDEV VARIANCE RANGE MINIMUM MAXIMUM MEAN MEDIAN MODE SKEWNESS SESKEW KURTOSIS SEKURT.

* Tạo multiple histograms theo nhóm
* (Giả sử có biến region: 1=Bắc, 2=Trung, 3=Nam)
GRAPH /HISTOGRAM=age 
  /PANEL ROWVAR=region ROWOP=CROSS.

// Tạo histogram trong Stata
// Tác giả: Thực hành kinh tế lượng
// Ví dụ: Phân tích phân phối tuổi lao động

* Xóa dữ liệu hiện tại
clear

* Nhập dữ liệu mẫu về tuổi
input age
36
25
38
46
55
68
72
55
36
38
67
45
22
48
91
46
52
61
58
55
end

* Tạo histogram cơ bản
histogram age, title("Phân phối tuổi lao động") ///
    xtitle("Tuổi") ytitle("Mật độ")

* Tạo histogram với overlay normal curve
histogram age, normal title("Phân phối tuổi với đường cong chuẩn") ///
    xtitle("Tuổi") ytitle("Mật độ")

* Tùy chỉnh số bins
histogram age, bin(6) title("Histogram với 6 bins") ///
    xtitle("Tuổi") ytitle("Mật độ")

* Histogram với density thay vì frequency
histogram age, density normal title("Density histogram") ///
    xtitle("Tuổi") ytitle("Mật độ")

* Thống kê mô tả
summarize age, detail

* Kiểm định tính chuẩn
swilk age
sktest age

* Tạo histogram theo nhóm (ví dụ theo region)
* generate region = ceil(_n/7)  // Tạo biến region giả định
* histogram age, by(region) title("Phân phối tuổi theo vùng")

* Lưu histogram
graph save "histogram_age.gph", replace
graph export "histogram_age.png", replace

# Tạo histogram trong R
# Tác giả: Thực hành kinh tế lượng
# Ví dụ: Phân tích phân phối tuổi lao động

# Cài đặt và load các package cần thiết
if (!require(ggplot2)) install.packages("ggplot2")
if (!require(dplyr)) install.packages("dplyr")
if (!require(nortest)) install.packages("nortest")
library(ggplot2)
library(dplyr)
library(nortest)

# Tạo dữ liệu mẫu về tuổi
age <- c(36, 25, 38, 46, 55, 68, 72, 55, 36, 38,
         67, 45, 22, 48, 91, 46, 52, 61, 58, 55)

# Tạo data frame
data <- data.frame(age = age)

# Histogram cơ bản với base R
hist(age, 
     main = "Phân phối tuổi lao động Việt Nam",
     xlab = "Tuổi", 
     ylab = "Tần số",
     col = "lightblue",
     border = "black")

# Histogram với đường cong normal
hist(age, 
     main = "Phân phối tuổi với đường cong chuẩn",
     xlab = "Tuổi", 
     ylab = "Tần số",
     col = "lightblue",
     border = "black",
     prob = TRUE)  # Hiển thị density thay vì frequency

# Thêm đường cong normal
curve(dnorm(x, mean = mean(age), sd = sd(age)), 
      add = TRUE, col = "red", lwd = 2)

# Histogram với ggplot2 (chuyên nghiệp hơn)
p1 <- ggplot(data, aes(x = age)) +
  geom_histogram(aes(y = ..density..), 
                 binwidth = 10, 
                 fill = "lightblue", 
                 color = "black", 
                 alpha = 0.7) +
  stat_function(fun = dnorm, 
                args = list(mean = mean(age), sd = sd(age)),
                color = "red", size = 1) +
  labs(title = "Phân phối tuổi lao động (ggplot2)",
       x = "Tuổi",
       y = "Mật độ") +
  theme_minimal()

print(p1)

# Histogram với số bins khác nhau
par(mfrow = c(2, 2))
hist(age, breaks = 5, main = "5 bins", col = "lightblue")
hist(age, breaks = 8, main = "8 bins", col = "lightgreen")
hist(age, breaks = 15, main = "15 bins", col = "pink")
hist(age, breaks = 20, main = "20 bins", col = "lightyellow")
par(mfrow = c(1, 1))

# Thống kê mô tả
summary(age)
cat("Độ lệch chuẩn:", sd(age), "\n")
cat("Skewness:", moments::skewness(age), "\n")
cat("Kurtosis:", moments::kurtosis(age), "\n")

# Kiểm định tính chuẩn
shapiro.test(age)
lillie.test(age)

# Tạo Q-Q plot để kiểm tra tính chuẩn
qqnorm(age, main = "Q-Q Plot")
qqline(age, col = "red")

# Histogram theo nhóm (ví dụ)
# Tạo biến region giả định
data$region <- rep(c("Bắc", "Trung", "Nam"), length.out = length(age))

ggplot(data, aes(x = age)) +
  geom_histogram(binwidth = 10, fill = "lightblue", alpha = 0.7) +
  facet_wrap(~region) +
  labs(title = "Phân phối tuổi theo vùng miền",
       x = "Tuổi",
       y = "Tần số") +
  theme_minimal()

# Lưu đồ thị
ggsave("histogram_age.png", width = 10, height = 6)

# Tạo histogram trong Python
# Tác giả: Thực hành kinh tế lượng
# Ví dụ: Phân tích phân phối tuổi lao động

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

# Thiết lập font tiếng Việt
plt.rcParams['font.family'] = 'DejaVu Sans'

# Tạo dữ liệu mẫu về tuổi
age_data = [36, 25, 38, 46, 55, 68, 72, 55, 36, 38,
           67, 45, 22, 48, 91, 46, 52, 61, 58, 55]

# Tạo DataFrame
df = pd.DataFrame({'age': age_data})

# Histogram cơ bản với matplotlib
plt.figure(figsize=(10, 6))
plt.hist(age_data, bins=8, color='lightblue', edgecolor='black', alpha=0.7)
plt.title('Phân phối tuổi lao động Việt Nam', fontsize=14, fontweight='bold')
plt.xlabel('Tuổi', fontsize=12)
plt.ylabel('Tần số', fontsize=12)
plt.grid(True, alpha=0.3)
plt.show()

# Histogram với đường cong normal
plt.figure(figsize=(10, 6))
count, bins, ignored = plt.hist(age_data, bins=8, density=True, 
                               color='lightblue', edgecolor='black', alpha=0.7)

# Thêm đường cong normal
mu, sigma = np.mean(age_data), np.std(age_data)
x = np.linspace(min(age_data), max(age_data), 100)
plt.plot(x, stats.norm.pdf(x, mu, sigma), 'r-', linewidth=2, 
         label=f'Normal(μ={mu:.1f}, σ={sigma:.1f})')

plt.title('Phân phối tuổi với đường cong chuẩn', fontsize=14, fontweight='bold')
plt.xlabel('Tuổi', fontsize=12)
plt.ylabel('Mật độ', fontsize=12)
plt.legend()
plt.grid(True, alpha=0.3)
plt.show()

# Histogram với seaborn (đẹp hơn)
plt.figure(figsize=(10, 6))
sns.histplot(data=df, x='age', bins=8, kde=True, stat='density')
plt.title('Phân phối tuổi với KDE (Seaborn)', fontsize=14, fontweight='bold')
plt.xlabel('Tuổi', fontsize=12)
plt.ylabel('Mật độ', fontsize=12)
plt.show()

# So sánh histogram với số bins khác nhau
fig, axes = plt.subplots(2, 2, figsize=(12, 10))
bins_list = [5, 8, 15, 20]
colors = ['lightblue', 'lightgreen', 'pink', 'lightyellow']

for i, (bins, color) in enumerate(zip(bins_list, colors)):
    ax = axes[i//2, i%2]
    ax.hist(age_data, bins=bins, color=color, edgecolor='black', alpha=0.7)
    ax.set_title(f'{bins} bins')
    ax.set_xlabel('Tuổi')
    ax.set_ylabel('Tần số')
    ax.grid(True, alpha=0.3)

plt.tight_layout()
plt.show()

# Thống kê mô tả
print("=== THỐNG KÊ MÔ TẢ ===")
print(f"Số quan sát: {len(age_data)}")
print(f"Trung bình: {np.mean(age_data):.2f}")
print(f"Trung vị: {np.median(age_data):.2f}")
print(f"Độ lệch chuẩn: {np.std(age_data, ddof=1):.2f}")
print(f"Min: {np.min(age_data)}")
print(f"Max: {np.max(age_data)}")
print(f"Skewness: {stats.skew(age_data):.3f}")
print(f"Kurtosis: {stats.kurtosis(age_data):.3f}")

# Kiểm định tính chuẩn
print("\n=== KIỂM ĐỊNH TÍNH CHUẨN ===")
# Shapiro-Wilk test
shapiro_stat, shapiro_p = shapiro(age_data)
print(f"Shapiro-Wilk test: W = {shapiro_stat:.4f}, p-value = {shapiro_p:.4f}")

# D'Agostino and Pearson's test
dagostino_stat, dagostino_p = normaltest(age_data)
print(f"D'Agostino test: χ² = {dagostino_stat:.4f}, p-value = {dagostino_p:.4f}")

if shapiro_p > 0.05:
    print("Kết luận: Dữ liệu có thể có phân phối chuẩn (p > 0.05)")
else:
    print("Kết luận: Dữ liệu không có phân phối chuẩn (p ≤ 0.05)")

# Q-Q plot
plt.figure(figsize=(8, 6))
stats.probplot(age_data, dist="norm", plot=plt)
plt.title('Q-Q Plot', fontsize=14, fontweight='bold')
plt.grid(True, alpha=0.3)
plt.show()

# Histogram theo nhóm (giả định có biến region)
np.random.seed(42)
regions = np.random.choice(['Bắc', 'Trung', 'Nam'], size=len(age_data))
df['region'] = regions

# Subplot theo region
fig, axes = plt.subplots(1, 3, figsize=(15, 5))
regions_unique = df['region'].unique()

for i, region in enumerate(regions_unique):
    data_region = df[df['region'] == region]['age']
    axes[i].hist(data_region, bins=5, color='lightblue', edgecolor='black', alpha=0.7)
    axes[i].set_title(f'Vùng {region}')
    axes[i].set_xlabel('Tuổi')
    axes[i].set_ylabel('Tần số')
    axes[i].grid(True, alpha=0.3)

plt.tight_layout()
plt.show()

# Boxplot để phát hiện outliers
plt.figure(figsize=(8, 6))
plt.boxplot(age_data, vert=True)
plt.title('Boxplot - Phát hiện outliers', fontsize=14, fontweight='bold')
plt.ylabel('Tuổi', fontsize=12)
plt.grid(True, alpha=0.3)
plt.show()

# Tính outliers bằng IQR method
Q1 = np.percentile(age_data, 25)
Q3 = np.percentile(age_data, 75)
IQR = Q3 - Q1
lower_bound = Q1 - 1.5 * IQR
upper_bound = Q3 + 1.5 * IQR

outliers = [x for x in age_data if x < lower_bound or x > upper_bound]
print(f"\n=== PHÁT HIỆN OUTLIERS ===")
print(f"Q1: {Q1}")
print(f"Q3: {Q3}")
print(f"IQR: {IQR}")
print(f"Lower bound: {lower_bound}")
print(f"Upper bound: {upper_bound}")
print(f"Outliers: {outliers}")

# Lưu đồ thị
plt.figure(figsize=(10, 6))
plt.hist(age_data, bins=8, color='lightblue', edgecolor='black', alpha=0.7)
plt.title('Phân phối tuổi lao động Việt Nam', fontsize=14, fontweight='bold')
plt.xlabel('Tuổi', fontsize=12)
plt.ylabel('Tần số', fontsize=12)
plt.grid(True, alpha=0.3)
plt.savefig('histogram_age.png', dpi=300, bbox_inches='tight')
plt.show()

print("\nĐồ thị đã được lưu thành file 'histogram_age.png'")