Analisis ini mengkaji pengaruh tingkat pendidikan terhadap
pendapatan individu menggunakan data cross-section
wage1 dari Wooldridge (2019). Dataset ini memuat informasi
upah per jam, tahun pendidikan, pengalaman kerja, dan karakteristik
individu lainnya untuk 526 pekerja di Amerika Serikat tahun
1976.
Metode yang digunakan adalah Ordinary Least Squares (OLS) dengan spesifikasi Mincer Earnings Equation — salah satu model paling fundamental dalam ekonometrika mikro untuk menganalisis hubungan antara modal manusia (human capital) dan pendapatan individu.
Ekonometrika Mikro adalah cabang ekonometrika yang berfokus pada analisis perilaku unit ekonomi individual — individu, rumah tangga, atau perusahaan. Metode utamanya menggunakan data cross-section, yaitu data yang dikumpulkan dari banyak individu pada satu titik waktu, untuk menjawab pertanyaan seperti:
Regresi OLS (Ordinary Least Squares) dengan data cross-section adalah metode inti dalam ekonometrika mikro. Metode ini termasuk dalam kerangka model regresi linier klasik, di mana estimator OLS bersifat BLUE (Best Linear Unbiased Estimator) apabila asumsi Gauss-Markov terpenuhi.
Studi ini berlandaskan pada Human Capital Theory yang dikembangkan oleh Gary Becker (1964) dan Jacob Mincer (1974). Teori ini menyatakan bahwa pendidikan merupakan bentuk investasi pada diri sendiri yang meningkatkan produktivitas, dan pada gilirannya meningkatkan upah di pasar tenaga kerja.
Mincer (1974) memformulasikan hubungan ini dalam persamaan yang kini dikenal sebagai Mincer Earnings Equation:
\[\ln(wage_i) = \beta_0 + \beta_1 \cdot educ_i + \beta_2 \cdot exper_i + \beta_3 \cdot exper_i^2 + \varepsilon_i\]
Penggunaan log upah sebagai variabel dependen memiliki dua alasan utama: (1) distribusi upah cenderung skewed ke kanan sehingga transformasi log menghasilkan distribusi yang lebih simetris, dan (2) koefisien \(\beta_1\) dapat langsung diinterpretasikan sebagai persentase kenaikan upah (return to education) untuk setiap tambahan satu tahun pendidikan.
Penambahan \(exper^2\) menangkap efek nonlinear pengalaman — upah meningkat seiring pengalaman tetapi pada tingkat yang semakin melambat (diminishing returns).
Regresi OLS cross-section dengan Mincer Equation tepat digunakan ketika:
library(ggplot2)
library(dplyr)
library(lmtest) # uji heteroskedastisitas (Breusch-Pagan)
library(car) # uji multikolinearitas (VIF)
library(nortest) # uji normalitas (Lilliefors)
library(sandwich) # robust standard errors
library(ggcorrplot) # korelasi plot
# Muat data wage1 dari file CSV (Wooldridge, 2019)
# Sumber: https://vincentarelbundock.github.io/Rdatasets/csv/wooldridge/wage1.csv
df <- read.csv("wage1.csv")
# Tampilkan struktur data
str(df)
## 'data.frame': 526 obs. of 25 variables:
## $ rownames: int 1 2 3 4 5 6 7 8 9 10 ...
## $ wage : num 3.1 3.24 3 6 5.3 ...
## $ educ : int 11 12 11 8 12 16 18 12 12 17 ...
## $ exper : int 2 22 2 44 7 9 15 5 26 22 ...
## $ tenure : int 0 2 0 28 2 8 7 3 4 21 ...
## $ nonwhite: int 0 0 0 0 0 0 0 0 0 0 ...
## $ female : int 1 1 0 0 0 0 0 1 1 0 ...
## $ married : int 0 1 0 1 1 1 0 0 0 1 ...
## $ numdep : int 2 3 2 0 1 0 0 0 2 0 ...
## $ smsa : int 1 1 0 1 0 1 1 1 1 1 ...
## $ northcen: int 0 0 0 0 0 0 0 0 0 0 ...
## $ south : int 0 0 0 0 0 0 0 0 0 0 ...
## $ west : int 1 1 1 1 1 1 1 1 1 1 ...
## $ construc: int 0 0 0 0 0 0 0 0 0 0 ...
## $ ndurman : int 0 0 0 0 0 0 0 0 0 0 ...
## $ trcommpu: int 0 0 0 0 0 0 0 0 0 0 ...
## $ trade : int 0 0 1 0 0 0 1 0 1 0 ...
## $ services: int 0 1 0 0 0 0 0 0 0 0 ...
## $ profserv: int 0 0 0 0 0 1 0 0 0 0 ...
## $ profocc : int 0 0 0 0 0 1 1 1 1 1 ...
## $ clerocc : int 0 0 0 1 0 0 0 0 0 0 ...
## $ servocc : int 0 1 0 0 0 0 0 0 0 0 ...
## $ lwage : num 1.13 1.18 1.1 1.79 1.67 ...
## $ expersq : int 4 484 4 1936 49 81 225 25 676 484 ...
## $ tenursq : int 0 4 0 784 4 64 49 9 16 441 ...
# Buat variabel baru
df <- df %>%
mutate(
lwage = log(wage), # log upah (variabel dependen)
exper2 = exper^2, # pengalaman kuadrat (nonlinearitas)
female = as.factor(female), # gender sebagai faktor
married = as.factor(married) # status menikah sebagai faktor
)
# Tampilkan beberapa baris pertama
head(df[, c("wage", "lwage", "educ", "exper", "exper2", "female", "married", "tenure")])
## wage lwage educ exper exper2 female married tenure
## 1 3.10 1.131402 11 2 4 1 0 0
## 2 3.24 1.175573 12 22 484 1 1 2
## 3 3.00 1.098612 11 2 4 0 0 0
## 4 6.00 1.791759 8 44 1936 0 1 28
## 5 5.30 1.667707 12 7 49 0 1 2
## 6 8.75 2.169054 16 9 81 0 1 8
# Statistik deskriptif variabel utama
df %>%
select(wage, lwage, educ, exper, tenure) %>%
summary()
## wage lwage educ exper
## Min. : 0.530 Min. :-0.6349 Min. : 0.00 Min. : 1.00
## 1st Qu.: 3.330 1st Qu.: 1.2030 1st Qu.:12.00 1st Qu.: 5.00
## Median : 4.650 Median : 1.5369 Median :12.00 Median :13.50
## Mean : 5.896 Mean : 1.6233 Mean :12.56 Mean :17.02
## 3rd Qu.: 6.880 3rd Qu.: 1.9286 3rd Qu.:14.00 3rd Qu.:26.00
## Max. :24.980 Max. : 3.2181 Max. :18.00 Max. :51.00
## tenure
## Min. : 0.000
## 1st Qu.: 0.000
## Median : 2.000
## Mean : 5.105
## 3rd Qu.: 7.000
## Max. :44.000
p1 <- ggplot(df, aes(x = wage)) +
geom_histogram(bins = 30, fill = "#2c7fb8", color = "white", alpha = 0.85) +
labs(title = "Distribusi Upah (Asli)", x = "Upah per Jam (USD)", y = "Frekuensi") +
theme_minimal(base_size = 12)
p2 <- ggplot(df, aes(x = lwage)) +
geom_histogram(bins = 30, fill = "#41b6c4", color = "white", alpha = 0.85) +
labs(title = "Distribusi Log(Upah)", x = "ln(Upah)", y = "Frekuensi") +
theme_minimal(base_size = 12)
gridExtra::grid.arrange(p1, p2, ncol = 2)
Distribusi upah: asli vs log-transformasi
Log-transformasi pada upah menghasilkan distribusi yang jauh lebih mendekati normal, sehingga lebih sesuai untuk estimasi OLS.
ggplot(df, aes(x = educ, y = lwage)) +
geom_jitter(alpha = 0.4, color = "#2c7fb8", width = 0.2) +
geom_smooth(method = "lm", color = "#e34a33", se = TRUE, linewidth = 1.2) +
labs(
title = "Pendidikan vs Log(Upah)",
subtitle = "Setiap titik mewakili satu individu",
x = "Tahun Pendidikan",
y = "ln(Upah per Jam)"
) +
theme_minimal(base_size = 12)
## `geom_smooth()` using formula = 'y ~ x'
Hubungan pendidikan dan log upah
model1 <- lm(lwage ~ educ, data = df)
summary(model1)
##
## Call:
## lm(formula = lwage ~ educ, data = df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.21158 -0.36393 -0.07263 0.29712 1.52339
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.583773 0.097336 5.998 3.74e-09 ***
## educ 0.082744 0.007567 10.935 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.4801 on 524 degrees of freedom
## Multiple R-squared: 0.1858, Adjusted R-squared: 0.1843
## F-statistic: 119.6 on 1 and 524 DF, p-value: < 2.2e-16
model2 <- lm(lwage ~ educ + exper + exper2, data = df)
summary(model2)
##
## Call:
## lm(formula = lwage ~ educ + exper + exper2, data = df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.96387 -0.29375 -0.04009 0.29497 1.30216
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.1279975 0.1059323 1.208 0.227
## educ 0.0903658 0.0074680 12.100 < 2e-16 ***
## exper 0.0410089 0.0051965 7.892 1.77e-14 ***
## exper2 -0.0007136 0.0001158 -6.164 1.42e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.4459 on 522 degrees of freedom
## Multiple R-squared: 0.3003, Adjusted R-squared: 0.2963
## F-statistic: 74.67 on 3 and 522 DF, p-value: < 2.2e-16
model3 <- lm(lwage ~ educ + exper + exper2 + female + married + tenure, data = df)
summary(model3)
##
## Call:
## lm(formula = lwage ~ educ + exper + exper2 + female + married +
## tenure, data = df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.81854 -0.25682 -0.02526 0.24755 1.18148
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.4158424 0.0993372 4.186 3.33e-05 ***
## educ 0.0798322 0.0068273 11.693 < 2e-16 ***
## exper 0.0300995 0.0051931 5.796 1.18e-08 ***
## exper2 -0.0006012 0.0001099 -5.472 6.95e-08 ***
## female1 -0.2911303 0.0362832 -8.024 6.88e-15 ***
## married1 0.0564494 0.0409259 1.379 0.168
## tenure 0.0160739 0.0028801 5.581 3.86e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.4014 on 519 degrees of freedom
## Multiple R-squared: 0.4361, Adjusted R-squared: 0.4296
## F-statistic: 66.91 on 6 and 519 DF, p-value: < 2.2e-16
# Perbandingan AIC dan R-squared ketiga model
cat("=== Perbandingan Model ===\n")
## === Perbandingan Model ===
cat(sprintf("Model 1 — R² = %.4f | AIC = %.2f\n", summary(model1)$r.squared, AIC(model1)))
## Model 1 — R² = 0.1858 | AIC = 724.76
cat(sprintf("Model 2 — R² = %.4f | AIC = %.2f\n", summary(model2)$r.squared, AIC(model2)))
## Model 2 — R² = 0.3003 | AIC = 649.06
cat(sprintf("Model 3 — R² = %.4f | AIC = %.2f\n", summary(model3)$r.squared, AIC(model3)))
## Model 3 — R² = 0.4361 | AIC = 541.51
cat("\n>> Model 3 dipilih sebagai model terbaik (R² tertinggi, AIC terendah)\n")
##
## >> Model 3 dipilih sebagai model terbaik (R² tertinggi, AIC terendah)
koef <- coef(model3)
cat("=== Interpretasi Koefisien Model 3 ===\n\n")
## === Interpretasi Koefisien Model 3 ===
cat(sprintf(
"1. PENDIDIKAN (educ): Setiap tambahan 1 tahun pendidikan meningkatkan upah sebesar %.1f%%\n",
koef["educ"] * 100
))
## 1. PENDIDIKAN (educ): Setiap tambahan 1 tahun pendidikan meningkatkan upah sebesar 8.0%
cat(sprintf(
"2. PENGALAMAN (exper): Pengalaman memiliki efek nonlinear (konkaf).\n Puncak upah dicapai pada pengalaman = %.1f tahun\n",
-koef["exper"] / (2 * koef["exper2"])
))
## 2. PENGALAMAN (exper): Pengalaman memiliki efek nonlinear (konkaf).
## Puncak upah dicapai pada pengalaman = 25.0 tahun
cat(sprintf(
"3. GENDER (female=1): Perempuan rata-rata memperoleh upah %.1f%% lebih rendah dari laki-laki\n",
(exp(koef["female1"]) - 1) * 100
))
## 3. GENDER (female=1): Perempuan rata-rata memperoleh upah -25.3% lebih rendah dari laki-laki
cat(sprintf(
"4. STATUS MENIKAH (married=1): Pekerja menikah memperoleh upah %.1f%% lebih tinggi\n",
(exp(koef["married1"]) - 1) * 100
))
## 4. STATUS MENIKAH (married=1): Pekerja menikah memperoleh upah 5.8% lebih tinggi
cat(sprintf(
"5. MASA KERJA (tenure): Setiap tambahan 1 tahun di perusahaan meningkatkan upah sebesar %.1f%%\n",
koef["tenure"] * 100
))
## 5. MASA KERJA (tenure): Setiap tambahan 1 tahun di perusahaan meningkatkan upah sebesar 1.6%
H₀: Residual berdistribusi normal
H₁: Residual tidak berdistribusi normal
resid_m3 <- residuals(model3)
# Plot QQ
par(mfrow = c(1, 2))
hist(resid_m3, breaks = 30, col = "#2c7fb8", border = "white",
main = "Histogram Residual", xlab = "Residual")
qqnorm(resid_m3, main = "Q-Q Plot Residual")
qqline(resid_m3, col = "#e34a33", lwd = 2)
par(mfrow = c(1, 1))
# Lilliefors test
lf <- lillie.test(resid_m3)
cat("Lilliefors Test — p-value:", round(lf$p.value, 4), "\n")
## Lilliefors Test — p-value: 0.01
if (lf$p.value > 0.05) {
cat(">> Residual berdistribusi normal. Asumsi normalitas TERPENUHI.\n")
} else {
cat(">> Residual tidak normal. Namun dengan n =", nrow(df),
"observasi, OLS tetap valid secara asimtotik (Central Limit Theorem).\n")
}
## >> Residual tidak normal. Namun dengan n = 526 observasi, OLS tetap valid secara asimtotik (Central Limit Theorem).
H₀: Variansi residual konstan (homoskedastis)
H₁: Variansi residual tidak konstan
(heteroskedastis)
# Plot residual vs fitted
plot(fitted(model3), resid_m3,
xlab = "Fitted Values", ylab = "Residual",
main = "Residual vs Fitted Values",
pch = 19, col = adjustcolor("#2c7fb8", alpha.f = 0.5))
abline(h = 0, col = "#e34a33", lty = 2, lwd = 2)
# Breusch-Pagan test
bp <- bptest(model3)
cat("Breusch-Pagan Test — p-value:", round(bp$p.value, 4), "\n")
## Breusch-Pagan Test — p-value: 0.0092
if (bp$p.value > 0.05) {
cat(">> Tidak ada heteroskedastisitas. Asumsi homoskedastisitas TERPENUHI.\n")
} else {
cat(">> Terdeteksi heteroskedastisitas. Digunakan Robust Standard Errors (HC3).\n")
cat("\n--- Hasil dengan Robust SE ---\n")
print(coeftest(model3, vcov = vcovHC(model3, type = "HC3")))
}
## >> Terdeteksi heteroskedastisitas. Digunakan Robust Standard Errors (HC3).
##
## --- Hasil dengan Robust SE ---
##
## t test of coefficients:
##
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.41584244 0.10835446 3.8378 0.0001394 ***
## educ 0.07983221 0.00780404 10.2296 < 2.2e-16 ***
## exper 0.03009952 0.00492995 6.1054 2.010e-09 ***
## exper2 -0.00060115 0.00010318 -5.8262 9.959e-09 ***
## female1 -0.29113027 0.03717593 -7.8312 2.743e-14 ***
## married1 0.05644940 0.04221616 1.3372 0.1817590
## tenure 0.01607388 0.00356680 4.5065 8.150e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Variance Inflation Factor (VIF) > 10 mengindikasikan masalah multikolinearitas serius.
vif_val <- vif(model3)
print(vif_val)
## educ exper exper2 female married tenure
## 1.164312 16.183797 14.924497 1.072362 1.302549 1.410458
cat("\n")
if (all(vif_val < 10)) {
cat(">> Semua VIF < 10. Tidak ada multikolinearitas serius. Asumsi TERPENUHI.\n")
} else {
cat(">> Ada variabel dengan VIF > 10. Perlu penanganan multikolinearitas.\n")
}
## >> Ada variabel dengan VIF > 10. Perlu penanganan multikolinearitas.
# Tidy koefisien untuk plotting
koef_df <- data.frame(
variabel = names(coef(model3))[-1],
estimate = coef(model3)[-1],
se = sqrt(diag(vcov(model3)))[-1]
) %>%
mutate(
lower = estimate - 1.96 * se,
upper = estimate + 1.96 * se,
signif = ifelse((lower > 0 & upper > 0) | (lower < 0 & upper < 0), "Signifikan", "Tidak Signifikan")
)
ggplot(koef_df, aes(x = reorder(variabel, estimate), y = estimate, color = signif)) +
geom_point(size = 3.5) +
geom_errorbar(aes(ymin = lower, ymax = upper), width = 0.2, linewidth = 1) +
geom_hline(yintercept = 0, linetype = "dashed", color = "gray50") +
scale_color_manual(values = c("Signifikan" = "#2c7fb8", "Tidak Signifikan" = "#bdbdbd")) +
coord_flip() +
labs(
title = "Koefisien Regresi dengan Confidence Interval 95%",
subtitle = "Model: ln(wage) ~ educ + exper + exper^2 + female + married + tenure",
x = NULL,
y = "Koefisien",
color = NULL
) +
theme_minimal(base_size = 12) +
theme(legend.position = "bottom")
Koefisien model dengan confidence interval 95%
# Prediksi untuk laki-laki menikah dengan pengalaman rata-rata
pred_df <- data.frame(
educ = seq(0, 18, by = 0.5),
exper = mean(df$exper),
exper2 = mean(df$exper)^2,
female = factor(0, levels = c(0, 1)),
married = factor(1, levels = c(0, 1)),
tenure = mean(df$tenure)
)
pred_df$lwage_pred <- predict(model3, newdata = pred_df)
ggplot(pred_df, aes(x = educ, y = lwage_pred)) +
geom_line(color = "#2c7fb8", linewidth = 1.5) +
geom_ribbon(aes(
ymin = predict(model3, newdata = pred_df, interval = "confidence")[, "lwr"],
ymax = predict(model3, newdata = pred_df, interval = "confidence")[, "upr"]
), alpha = 0.15, fill = "#2c7fb8") +
labs(
title = "Return to Education: Prediksi Log Upah vs Tahun Pendidikan",
subtitle = "Holding constant: pengalaman & tenure rata-rata, laki-laki, menikah",
x = "Tahun Pendidikan",
y = "Prediksi ln(Upah per Jam)"
) +
theme_minimal(base_size = 12)
Prediksi log upah berdasarkan tahun pendidikan
Analisis regresi OLS terhadap data wage1 (Wooldridge,
2019) menggunakan spesifikasi Mincer Earnings Equation menghasilkan
beberapa temuan utama:
Return to Education — Setiap tambahan satu tahun pendidikan meningkatkan upah sebesar ±8%, konsisten dengan Human Capital Theory. Hasil ini signifikan secara statistik pada α = 1%.
Efek Pengalaman Nonlinear — Pengalaman kerja memiliki efek positif yang menurun (diminishing returns), dengan puncak upah dicapai sekitar 25 tahun pengalaman.
Gender Wage Gap — Perempuan memperoleh upah rata-rata ±25.3% lebih rendah dari laki-laki dengan karakteristik serupa, mengindikasikan adanya diskriminasi upah berbasis gender.
Uji Asumsi — Model memenuhi sebagian besar asumsi klasik OLS. Apabila terdeteksi heteroskedastisitas, Robust Standard Errors (HC3) digunakan untuk koreksi inferensi.
Kelayakan Model — Model 3 (Mincer + kontrol demografis) memiliki R² tertinggi dan AIC terendah dibanding spesifikasi yang lebih sederhana, menjadikannya model terpilih.
Referensi:
wage1 — tersedia di https://vincentarelbundock.github.io/Rdatasets/csv/wooldridge/wage1.csv