cat("第1列：全样本真实婚姻匹配\n")
cat("年龄相关系数: ", cor(real_husband$husband_age, real_wife$wife_age, use = "complete.obs"), "\n")
cat("收入相关系数: ", cor(real_husband$husband_inc, real_wife$wife_inc, use = "complete.obs"), "\n")
cat("教育相关系数: ", cor(real_husband$husband_ysh, real_wife$wife_ysh, use = "complete.obs"), "\n\n")
cat("第2列：双收入家庭真实婚姻匹配\n")
cat("年龄相关系数: ", cor(married_sample_di$rage, married_sample_di$sage, use = "complete.obs"), "\n")
cat("收入相关系数: ", cor(married_sample_di$rinc, married_sample_di$sinc, use = "complete.obs"), "\n")
cat("教育相关系数: ", cor(married_sample_di$rysh, married_sample_di$sysh, use = "complete.obs"), "\n\n")
cat("第3列：全样本男性最优匹配\n")
cat("年龄：", mean(table3_man[[1]]), " (SD=", sd(table3_man[[1]]), ")\n")
cat("收入：", mean(table3_man[[2]]), " (SD=", sd(table3_man[[2]]), ")\n")
cat("教育：", mean(table3_man[[3]]), " (SD=", sd(table3_man[[3]]), ")\n\n")
cat("第4列：双收入家庭男性最优匹配\n")
cat("年龄：", mean(tableA21_man_di[[1]]), " (SD=", sd(tableA21_man_di[[1]]), ")\n")
cat("收入：", mean(tableA21_man_di[[2]]), " (SD=", sd(tableA21_man_di[[2]]), ")\n")
cat("教育：", mean(tableA21_man_di[[3]]), " (SD=", sd(tableA21_man_di[[3]]), ")\n\n")
cat("第5列：全样本女性最优匹配\n")
cat("年龄：", mean(table3_woman[[1]]), " (SD=", sd(table3_woman[[1]]), ")\n")
cat("收入：", mean(table3_woman[[2]]), " (SD=", sd(table3_woman[[2]]), ")\n")
cat("教育：", mean(table3_woman[[3]]), " (SD=", sd(table3_woman[[3]]), ")\n\n")
cat("第6列：双收入家庭女性最优匹配\n")
cat("年龄：", mean(tableA21_woman_di[[1]]), " (SD=", sd(tableA21_woman_di[[1]]), ")\n")
cat("收入：", mean(tableA21_woman_di[[2]]), " (SD=", sd(tableA21_woman_di[[2]]), ")\n")
cat("教育：", mean(tableA21_woman_di[[3]]), " (SD=", sd(tableA21_woman_di[[3]]), ")\n\n")
sink()
#创建真实丈夫的样本
real_husband = data.frame(matrix(nrow =nrow(married_sample)))
real_husband$Mid = married_sample$id
real_husband$husband_age = ifelse(married_sample$male=="男性",married_sample$rage,
married_sample$sage)
real_husband$husband_inc = ifelse(married_sample$male=="男性",married_sample$rinc,
married_sample$sinc)
real_husband$husband_ysh = ifelse(married_sample$male=="男性",married_sample$rysh,
married_sample$sysh)
real_husband$Marstatus = ifelse(married_sample$male=="男性",married_sample$MarStatus,
married_sample$MarStatus)
real_husband = unique(real_husband)
# 如果married_sample中male=="男性"的id在real_wife数据框中id列里有，则生成一列response==1，否则为0
real_husband$response <- ifelse(real_husband$Mid %in% married_sample$id[married_sample$male == "男性"], 1, 0)
#real_husband生成新id,受访者保持不变
real_husband$id = ifelse(real_husband$response==1,real_husband$Mid,real_husband$Mid+8148)
#生成id
HW_pref$id[10:2237] = as.numeric(gsub("[^0-9]", "", rownames(HW_pref)[10:2237]))
#匹配id,这样每个受访者男性的固定效果和std均已生成
real_husband = merge(real_husband,HW_pref[10:2237,],by.x = "id",by.y = "id",all.x = T)
#固定效果回归男性们的年龄、教育、收入,强制过原点
reservation_HW = list(summary(lm(HW_coef~0+husband_age+I(husband_inc/10000)+husband_ysh,data = real_husband[real_husband$response==1,])),
summary(lm(HW_std ~0+husband_age+I(husband_inc/10000)+husband_ysh,data = real_husband[real_husband$response==1,])))
real_husband$HW_coef[real_husband$response==0] =
reservation_HW[[1]][["coefficients"]][1,1]*real_husband$husband_age+
reservation_HW[[1]][["coefficients"]][2,1]*real_husband$husband_inc/10000+
reservation_HW[[1]][["coefficients"]][3,1]*real_husband$husband_ysh
real_husband$HW_std[real_husband$response==0] =
reservation_HW[[2]][["coefficients"]][1,1]*real_husband$husband_age+
reservation_HW[[2]][["coefficients"]][2,1]*real_husband$husband_inc/10000+
reservation_HW[[2]][["coefficients"]][3,1]*real_husband$husband_ysh
#创建真实妻子的样本
real_wife = data.frame(matrix(nrow =nrow(married_sample)))
real_wife$Fid = married_sample$id
real_wife$wife_age = ifelse(married_sample$male=="女性",married_sample$rage,
married_sample$sage)
real_wife$wife_inc = ifelse(married_sample$male=="女性",married_sample$rinc,
married_sample$sinc)
real_wife$wife_ysh = ifelse(married_sample$male=="女性",married_sample$rysh,
married_sample$sysh)
real_wife$Marstatus = ifelse(married_sample$male=="女性",married_sample$MarStatus,
married_sample$MarStatus)
real_wife = unique(real_wife)
# 如果married_sample中male=="女性"的id在real_wife数据框中id列里有，则生成一列response==1，否则为0
real_wife$response <- ifelse(real_wife$Fid %in% married_sample$id[married_sample$male == "女性"], 1, 0)
#real_wife生成新id,受访者保持不变
real_wife$id = ifelse(real_wife$response==1,real_wife$Fid,real_wife$Fid+8148)
#生成id
WH_pref$id[10:2745] = as.numeric(gsub("[^0-9]", "", rownames(WH_pref)[10:2745]))
#匹配id,这样每个受访者女性的固定效果和std均已生成
real_wife = merge(real_wife,WH_pref[10:2745,],by.x = "id",by.y = "id",all.x = T)
#固定效果回归女性们的年龄、教育、收入,强制过原点
reservation_WH = list(summary(lm(WH_coef~0+wife_age+I(wife_inc/10000)+wife_ysh,data = real_wife[real_wife$response==1,])),
summary(lm(WH_std ~0+wife_age+I(wife_inc/10000)+wife_ysh,data = real_wife[real_wife$response==1,])))
real_wife$WH_coef[real_wife$response==0] =
reservation_WH[[1]][["coefficients"]][1,1]*real_wife$wife_age+
reservation_WH[[1]][["coefficients"]][2,1]*real_wife$wife_inc/10000+
reservation_WH[[1]][["coefficients"]][3,1]*real_wife$wife_ysh
real_wife$WH_std[real_wife$response==0] =
reservation_WH[[2]][["coefficients"]][1,1]*real_wife$wife_age+
reservation_WH[[2]][["coefficients"]][2,1]*real_wife$wife_inc/10000+
reservation_WH[[2]][["coefficients"]][3,1]*real_wife$wife_ysh
#Mid 和 Fid 一样的是一家子
real_husband=real_husband[order(real_husband$Mid),]
real_wife=real_wife[order(real_wife$Fid),]
set.seed(7)
random_id = sample(real_husband$Mid,500)
real_husband = real_husband[real_husband$Mid%in%random_id,]
real_wife = real_wife[real_wife$Fid%in%random_id,]
#按照真实夫妻顺序依次排列
real_husband=real_husband[order(real_husband$Mid),]
real_wife=real_wife[order(real_wife$Fid),]
#dif_MF是个list, 第一个元素是500男和1号女的年龄差等变量，第二个元素是500男和2号女的年龄差等。
dif_MF_n = data.frame(matrix(nrow=nrow(real_husband),ncol=6))
#list总共500个
dif_MF = list()
#dif_MF是所有男性和第j个女性的关系
for(j in 1:nrow(real_wife)){
dif_MF[[j]] = dif_MF_n
}
for (j in 1:nrow(real_wife)) {
if(j %%50==0)
print(paste("Iteration",j))
for (n in 1:nrow(real_husband)) {
dif_MF[[j]][n,1] = ifelse(real_husband$husband_age[n]>real_wife$wife_age[j],
real_husband$husband_age[n]-real_wife$wife_age[j],0)
dif_MF[[j]][n,2] = ifelse(real_husband$husband_age[n]<real_wife$wife_age[j],
real_wife$wife_age[j]-real_husband$husband_age[n],0)
dif_MF[[j]][n,3] = ifelse(real_husband$husband_inc[n]>real_wife$wife_inc[j],
(real_husband$husband_inc[n]+1)/(real_wife$wife_inc[j]+1)-1,0)
dif_MF[[j]][n,4] = ifelse(real_husband$husband_inc[n]<real_wife$wife_inc[j],
1-(real_husband$husband_inc[n]+1)/(real_wife$wife_inc[j]+1),0)
dif_MF[[j]][n,5] = ifelse(real_husband$husband_ysh[n]>real_wife$wife_ysh[j],
real_husband$husband_ysh[n]-real_wife$wife_ysh[j],0)
dif_MF[[j]][n,6] = ifelse(real_husband$husband_ysh[n]<real_wife$wife_ysh[j],
real_wife$wife_ysh[j]-real_husband$husband_ysh[n],0)
}
}
###############################################################
#dif_FM是个list, 第一个元素是100女和1号男的年龄差等变量，第二个元素是100女和2号男的年龄差等。
dif_FM_n = data.frame(matrix(nrow=nrow(real_wife),ncol=6))
dif_FM = list()
for(n in 1:nrow(real_husband)){
dif_FM[[n]] = dif_FM_n
}
for (j in 1:nrow(real_husband)) {
if(j %%50==0)
print(paste("Iteration",j))
for (n in 1:nrow(real_wife)) {
dif_FM[[j]][n,1] = ifelse(real_wife$wife_age[n]>real_husband$husband_age[j],
real_wife$wife_age[n]-real_husband$husband_age[j],0)
dif_FM[[j]][n,2] = ifelse(real_wife$wife_age[n]<real_husband$husband_age[j],
-(real_wife$wife_age[n]-real_husband$husband_age[j]),0)
dif_FM[[j]][n,3] = ifelse(real_wife$wife_inc[n]>real_husband$husband_inc[j],
(real_wife$wife_inc[n]+1)/(real_husband$husband_inc[j]+1)-1,0)
dif_FM[[j]][n,4] = ifelse(real_wife$wife_inc[n]<real_husband$husband_inc[j],
1-(real_wife$wife_inc[n]+1)/(real_husband$husband_inc[j]+1),0)
dif_FM[[j]][n,5] = ifelse(real_wife$wife_ysh[n]>real_husband$husband_ysh[j],
real_wife$wife_ysh[n]-real_husband$husband_ysh[j],0)
dif_FM[[j]][n,6] = ifelse(real_wife$wife_ysh[n]<real_husband$husband_ysh[j],
-(real_wife$wife_ysh[n]-real_husband$husband_ysh[j]),0)
}
}
HWcom_pref = readRDS("HWcom_pref.rds")
WHcom_pref = readRDS("WHcom_pref.rds")
View(WHcom_pref)
WHcom_pref[[1]]
#让偏好11列为0，代表随机效用项=0
WHcom_pref = lapply(WHcom_pref, function(mat) {
mat[, c(11)]=0
return(mat)})
HWcom_pref = lapply(HWcom_pref, function(mat) {
mat[, c(11)]=0
return(mat)})
WHcom_pref[[1]]
#两个矩阵相乘，得出每位女性给每位男性的打分,dif_MF[[n]]是第n号女性的每位男生比她大的年龄等矩阵。
#WHcom_pref[[r]][n]是玩第r轮游戏，第n个女性的所有参数  （beta1,beta2,,,,,alphai,ui）'
#前面的矩阵长这样“（age男1,age+男1-女1，age-男1-女1，，，，1，1）第二行是第二个男性和第一个女性的关系
WH_all_r = matrix(ncol = nrow(real_husband),nrow=nrow(real_wife))
WH_all = list()
#WH_all是每位女性对所有潜在配偶的排名，100轮游戏的结果，第一个元素是个5741*5741的矩阵，第一列是第一位女性对所有男性的打分
for (round in 1:100) WH_all[[round]] = WH_all_r
for(round in 1:100){
#第一列是第一位女性对所有男人的效用值，第二列是第二女性
if(round %%10==0)
print(paste("Iteration",round))
for (n in 1:nrow(real_wife)){
WH_all[[round]][,n] = as.matrix(cbind(real_husband$husband_age,(dif_MF[[n]][,1:2])^2,
real_husband$husband_inc/10000,(dif_MF[[n]][,3:4])^2,
real_husband$husband_ysh,(dif_MF[[n]][,5:6])^2,1,1)) %*% matrix(WHcom_pref[[round]][n,])
}
}
sort_mat <- function(mat) {
mat_sorted <- sapply(seq_len(ncol(mat)), function(col_idx) {
sorted_col <- c(order(-mat[, col_idx],na.last = NA),
rep(NA,nrow(real_husband)-length(order(-mat[, col_idx],na.last = NA))))
})
return(mat_sorted)
}
#WH_all_pref[[1]],第一轮游戏中每位女性对所有男性排序
WH_all_pref = lapply(WH_all, sort_mat)
####################################################################################################
HW_all_r = matrix(ncol = nrow(real_wife),nrow=nrow(real_husband))
HW_all = list()
#HW_all是每位男性对所有潜在配偶的排名，100轮游戏的结果，第一个元素是个100*100的矩阵，第一列是第一位男性对所有女性的打分
for (round in 1:100) HW_all[[round]] = HW_all_r
for(round in 1:100){
#第一列是第一位女性对所有男人的效用值，第二列是第二女性
if(round %%10==0)
print(paste("Iteration",round))
for (n in 1:nrow(real_husband)){
HW_all[[round]][,n] = as.matrix(cbind(real_wife$wife_age,(dif_FM[[n]][,1:2])^2,
real_wife$wife_inc/10000,(dif_FM[[n]][,3:4])^2,
real_wife$wife_ysh,(dif_FM[[n]][,5:6])^2,1,1)) %*% matrix(HWcom_pref[[round]][n,])
}
}
#HW_all_pref[[1]],第一轮游戏中每位女性对所有男性排序
HW_all_pref = lapply(HW_all, sort_mat)
#求GS算法
corr_age = c()
corr_inc = c()
corr_ysh = c()
diff_age = c()
diff_inc = c()
diff_ysh = c()
std_age = c()
std_inc = c()
std_ysh = c()
nmatch =c()
for (round in 1:100) {
if(round %%10==0)
print(paste("Iteration",round))
#s给男性的偏好，算man-optimal
s.prefs <- HW_all_pref[[round]]
c.prefs <- WH_all_pref[[round]]
matching = iaa(s.prefs = s.prefs, c.prefs = c.prefs,
acceptance="deferred")$matchings
corr_age[round] = cor(real_husband$husband_age[matching$college],real_wife$wife_age[matching$student],use = "complete.obs")
corr_inc[round] = cor(real_husband$husband_inc[matching$college],real_wife$wife_inc[matching$student],use = "complete.obs")
corr_ysh[round] = cor(real_husband$husband_ysh[matching$college],real_wife$wife_ysh[matching$student],use = "complete.obs")
diff_age[round] = mean(real_husband$husband_age[matching$college] - real_wife$wife_age[matching$student],na.rm = T)
diff_inc[round] = mean(real_husband$husband_inc[matching$college] - real_wife$wife_inc[matching$student],na.rm = T)
diff_ysh[round] = mean(real_husband$husband_ysh[matching$college] - real_wife$wife_ysh[matching$student],na.rm = T)
std_age[round] = sd(real_husband$husband_age[matching$college] - real_wife$wife_age[matching$student],na.rm = T)
std_inc[round] = sd(real_husband$husband_inc[matching$college] - real_wife$wife_inc[matching$student],na.rm = T)
std_ysh[round] = sd(real_husband$husband_ysh[matching$college] - real_wife$wife_ysh[matching$student],na.rm = T)
nmatch[round] = nrow(matching)
}
mean(corr_age)
sd(corr_age)
mean(corr_inc)
sd(corr_inc)
mean(corr_ysh)
sd(corr_ysh)
mean(diff_age)
mean(diff_inc)
mean(diff_ysh)
mean(std_age)
mean(std_inc)
mean(std_ysh)
tableA22_man=list(corr_age,corr_inc,corr_ysh,diff_age,diff_inc,diff_ysh,std_age,std_inc,std_ysh)
WHcom_pref = readRDS("WHcom_pref.rds")
WHcom_pref_r = matrix(nrow=nrow(real_wife))
WHcom_pref = list()
#生成100次的女性偏好矩阵
for(round in 1:100){
WHcom_pref[[round]] = WHcom_pref_r
}
#给每个矩阵里填数字，每个矩阵的1：9列是slope,最后一列是error
for (round in 1:100) {
WHcom_pref[[round]] =  cbind(mvrnorm(nrow(real_wife),married_WH[["coefficients"]][1:9,1],vcov(married_WH)[1:9,1:9]),
real_wife$WH_coef,
rnorm(nrow(real_wife),0,sd(resid(married_WH))))
}
#让偏好11列为0，代表随机效用项=0
WHcom_pref = lapply(WHcom_pref, function(mat) {
mat[, c(11)]=0
return(mat)})
HWcom_pref = readRDS("HWcom_pref.rds")
HWcom_pref = lapply(HWcom_pref, function(mat) {
mat[, c(11)]=0
return(mat)})
WH_all_r = matrix(ncol = nrow(real_husband),nrow=nrow(real_wife))
WH_all = list()
#WH_all是每位女性对所有潜在配偶的排名，100轮游戏的结果，第一个元素是个5741*5741的矩阵，第一列是第一位女性对所有男性的打分
for (round in 1:100) WH_all[[round]] = WH_all_r
for(round in 1:100){
#第一列是第一位女性对所有男人的效用值，第二列是第二女性
if(round %%10==0)
print(paste("Iteration",round))
for (n in 1:nrow(real_wife)){
WH_all[[round]][,n] = as.matrix(cbind(real_husband$husband_age,(dif_MF[[n]][,1:2])^2,
real_husband$husband_inc/10000,(dif_MF[[n]][,3:4])^2,
real_husband$husband_ysh,(dif_MF[[n]][,5:6])^2,1,1)) %*% matrix(WHcom_pref[[round]][n,])
}
}
sort_mat <- function(mat) {
mat_sorted <- sapply(seq_len(ncol(mat)), function(col_idx) {
sorted_col <- c(order(-mat[, col_idx],na.last = NA),
rep(NA,nrow(real_husband)-length(order(-mat[, col_idx],na.last = NA))))
})
return(mat_sorted)
}
#WH_all_pref[[1]],第一轮游戏中每位女性对所有男性排序
WH_all_pref = lapply(WH_all, sort_mat)
####################################################################################################
HW_all_r = matrix(ncol = nrow(real_wife),nrow=nrow(real_husband))
HW_all = list()
#HW_all是每位男性对所有潜在配偶的排名，100轮游戏的结果，第一个元素是个100*100的矩阵，第一列是第一位男性对所有女性的打分
for (round in 1:100) HW_all[[round]] = HW_all_r
for(round in 1:100){
#第一列是第一位女性对所有男人的效用值，第二列是第二女性
if(round %%10==0)
print(paste("Iteration",round))
for (n in 1:nrow(real_husband)){
HW_all[[round]][,n] = as.matrix(cbind(real_wife$wife_age,(dif_FM[[n]][,1:2])^2,
real_wife$wife_inc/10000,(dif_FM[[n]][,3:4])^2,
real_wife$wife_ysh,(dif_FM[[n]][,5:6])^2,1,1)) %*% matrix(HWcom_pref[[round]][n,])
}
}
#HW_all_pref[[1]],第一轮游戏中每位女性对所有男性排序
HW_all_pref = lapply(HW_all, sort_mat)
#求GS算法
corr_age = c()
corr_inc = c()
corr_ysh = c()
diff_age = c()
diff_inc = c()
diff_ysh = c()
std_age = c()
std_inc = c()
std_ysh = c()
nmatch =c()
for (round in 1:100) {
if(round %%10==0)
print(paste("Iteration",round))
#s给男性的偏好，算man-optimal
s.prefs <- HW_all_pref[[round]]
c.prefs <- WH_all_pref[[round]]
matching = iaa(s.prefs = s.prefs, c.prefs = c.prefs,
acceptance="deferred")$matchings
corr_age[round] = cor(real_husband$husband_age[matching$college],real_wife$wife_age[matching$student],use = "complete.obs")
corr_inc[round] = cor(real_husband$husband_inc[matching$college],real_wife$wife_inc[matching$student],use = "complete.obs")
corr_ysh[round] = cor(real_husband$husband_ysh[matching$college],real_wife$wife_ysh[matching$student],use = "complete.obs")
diff_age[round] = mean(real_husband$husband_age[matching$college] - real_wife$wife_age[matching$student],na.rm = T)
diff_inc[round] = mean(real_husband$husband_inc[matching$college] - real_wife$wife_inc[matching$student],na.rm = T)
diff_ysh[round] = mean(real_husband$husband_ysh[matching$college] - real_wife$wife_ysh[matching$student],na.rm = T)
std_age[round] = sd(real_husband$husband_age[matching$college] - real_wife$wife_age[matching$student],na.rm = T)
std_inc[round] = sd(real_husband$husband_inc[matching$college] - real_wife$wife_inc[matching$student],na.rm = T)
std_ysh[round] = sd(real_husband$husband_ysh[matching$college] - real_wife$wife_ysh[matching$student],na.rm = T)
nmatch[round] = nrow(matching)
}
mean(corr_age)
sd(corr_age)
mean(corr_inc)
sd(corr_inc)
mean(corr_ysh)
sd(corr_ysh)
mean(diff_age)
mean(diff_inc)
mean(diff_ysh)
mean(std_age)
mean(std_inc)
mean(std_ysh)
tableA22_man=list(corr_age,corr_inc,corr_ysh,diff_age,diff_inc,diff_ysh,std_age,std_inc,std_ysh)
#求GS算法
corr_age = c()
corr_inc = c()
corr_ysh = c()
diff_age = c()
diff_inc = c()
diff_ysh = c()
std_age = c()
std_inc = c()
std_ysh = c()
nmatch =c()
for (round in 1:100) {
if(round %%10==0)
print(paste("Iteration",round))
#s给男性的偏好，算man-optimal
s.prefs <- WH_all_pref[[round]]
c.prefs <- HW_all_pref[[round]]
matching = iaa(s.prefs = s.prefs, c.prefs = c.prefs,
acceptance="deferred")$matchings
corr_age[round] = cor(real_husband$husband_age[matching$college],real_wife$wife_age[matching$student],use = "complete.obs")
corr_inc[round] = cor(real_husband$husband_inc[matching$college],real_wife$wife_inc[matching$student],use = "complete.obs")
corr_ysh[round] = cor(real_husband$husband_ysh[matching$college],real_wife$wife_ysh[matching$student],use = "complete.obs")
diff_age[round] = mean(real_husband$husband_age[matching$college] - real_wife$wife_age[matching$student],na.rm = T)
diff_inc[round] = mean(real_husband$husband_inc[matching$college] - real_wife$wife_inc[matching$student],na.rm = T)
diff_ysh[round] = mean(real_husband$husband_ysh[matching$college] - real_wife$wife_ysh[matching$student],na.rm = T)
std_age[round] = sd(real_husband$husband_age[matching$college] - real_wife$wife_age[matching$student],na.rm = T)
std_inc[round] = sd(real_husband$husband_inc[matching$college] - real_wife$wife_inc[matching$student],na.rm = T)
std_ysh[round] = sd(real_husband$husband_ysh[matching$college] - real_wife$wife_ysh[matching$student],na.rm = T)
nmatch[round] = nrow(matching)
}
mean(corr_age)
sd(corr_age)
mean(corr_inc)
sd(corr_inc)
mean(corr_ysh)
sd(corr_ysh)
mean(diff_age)
mean(diff_inc)
mean(diff_ysh)
mean(std_age)
mean(std_inc)
mean(std_ysh)
tableA22_woman=list(corr_age,corr_inc,corr_ysh,diff_age,diff_inc,diff_ysh,std_age,std_inc,std_ysh)
sink("附录Table2-2_log.txt")
cat("附录表2-2 匹配稳定性对比\n\n")
cat("第1列：男性最优\n")
cat("年龄：", mean(table3_man[[1]]), " (SD=", sd(table3_man[[1]]), ")\n")
cat("收入：", mean(table3_man[[2]]), " (SD=", sd(table3_man[[2]]), ")\n")
cat("教育：", mean(table3_man[[3]]), " (SD=", sd(table3_man[[3]]), ")\n\n")
cat("第2列：男性最优（随机效用=0）\n")
cat("年龄：", mean(tableA22_man[[1]]), " (SD=", sd(tableA22_man[[1]]), ")\n")
cat("收入：", mean(tableA22_man[[2]]), " (SD=", sd(tableA22_man[[2]]), ")\n")
cat("教育：", mean(tableA22_man[[3]]), " (SD=", sd(tableA22_man[[3]]), ")\n\n")
cat("第3列：女性最优\n")
cat("年龄：", mean(table3_woman[[1]]), " (SD=", sd(table3_woman[[1]]), ")\n")
cat("收入：", mean(table3_woman[[2]]), " (SD=", sd(table3_woman[[2]]), ")\n")
cat("教育：", mean(table3_woman[[3]]), " (SD=", sd(table3_woman[[3]]), ")\n\n")
cat("第4列：女性最优（随机效用=0）\n")
cat("年龄：", mean(tableA22_woman[[1]]), " (SD=", sd(tableA22_woman[[1]]), ")\n")
cat("收入：", mean(tableA22_woman[[2]]), " (SD=", sd(tableA22_woman[[2]]), ")\n")
cat("教育：", mean(tableA22_woman[[3]]), " (SD=", sd(tableA22_woman[[3]]), ")\n\n")
sink()
save.image()
library(readstata13)
library(tidyverse)
library(ggpubr)
library(MASS)
library(matchingMarkets)
library(bife)
#创建真实丈夫的样本
real_husband = data.frame(matrix(nrow =nrow(married_sample)))
real_husband$Mid = married_sample$id
real_husband$husband_age = ifelse(married_sample$male=="男性",married_sample$rage,
married_sample$sage)
real_husband$husband_inc = ifelse(married_sample$male=="男性",married_sample$rinc,
married_sample$sinc)
real_husband$husband_ysh = ifelse(married_sample$male=="男性",married_sample$rysh,
married_sample$sysh)
real_husband$Marstatus = ifelse(married_sample$male=="男性",married_sample$MarStatus,
married_sample$MarStatus)
real_husband = unique(real_husband)
# 如果married_sample中male=="男性"的id在real_wife数据框中id列里有，则生成一列response==1，否则为0
real_husband$response <- ifelse(real_husband$Mid %in% married_sample$id[married_sample$male == "男性"], 1, 0)
#real_husband生成新id,受访者保持不变
real_husband$id = ifelse(real_husband$response==1,real_husband$Mid,real_husband$Mid+8148)
#生成id
HW_pref$id[10:2237] = as.numeric(gsub("[^0-9]", "", rownames(HW_pref)[10:2237]))
#匹配id,这样每个受访者男性的固定效果和std均已生成
real_husband = merge(real_husband,HW_pref[10:2237,],by.x = "id",by.y = "id",all.x = T)
#固定效果回归男性们的年龄、教育、收入,强制过原点
reservation_HW = list(summary(lm(HW_coef~0+husband_age+I(husband_inc/10000)+husband_ysh,data = real_husband[real_husband$response==1,])),
summary(lm(HW_std ~0+husband_age+I(husband_inc/10000)+husband_ysh,data = real_husband[real_husband$response==1,])))
real_husband$HW_coef[real_husband$response==0] =
reservation_HW[[1]][["coefficients"]][1,1]*real_husband$husband_age+
reservation_HW[[1]][["coefficients"]][2,1]*real_husband$husband_inc/10000+
reservation_HW[[1]][["coefficients"]][3,1]*real_husband$husband_ysh
real_husband$HW_std[real_husband$response==0] =
reservation_HW[[2]][["coefficients"]][1,1]*real_husband$husband_age+
reservation_HW[[2]][["coefficients"]][2,1]*real_husband$husband_inc/10000+
reservation_HW[[2]][["coefficients"]][3,1]*real_husband$husband_ysh
#创建真实妻子的样本
real_wife = data.frame(matrix(nrow =nrow(married_sample)))
real_wife$Fid = married_sample$id
real_wife$wife_age = ifelse(married_sample$male=="女性",married_sample$rage,
married_sample$sage)
real_wife$wife_inc = ifelse(married_sample$male=="女性",married_sample$rinc,
married_sample$sinc)
real_wife$wife_ysh = ifelse(married_sample$male=="女性",married_sample$rysh,
married_sample$sysh)
real_wife$Marstatus = ifelse(married_sample$male=="女性",married_sample$MarStatus,
married_sample$MarStatus)
real_wife = unique(real_wife)
# 如果married_sample中male=="女性"的id在real_wife数据框中id列里有，则生成一列response==1，否则为0
real_wife$response <- ifelse(real_wife$Fid %in% married_sample$id[married_sample$male == "女性"], 1, 0)
#real_wife生成新id,受访者保持不变
real_wife$id = ifelse(real_wife$response==1,real_wife$Fid,real_wife$Fid+8148)
#生成id
WH_pref$id[10:2745] = as.numeric(gsub("[^0-9]", "", rownames(WH_pref)[10:2745]))
#匹配id,这样每个受访者女性的固定效果和std均已生成
real_wife = merge(real_wife,WH_pref[10:2745,],by.x = "id",by.y = "id",all.x = T)
#固定效果回归女性们的年龄、教育、收入,强制过原点
reservation_WH = list(summary(lm(WH_coef~0+wife_age+I(wife_inc/10000)+wife_ysh,data = real_wife[real_wife$response==1,])),
summary(lm(WH_std ~0+wife_age+I(wife_inc/10000)+wife_ysh,data = real_wife[real_wife$response==1,])))
real_wife$WH_coef[real_wife$response==0] =
reservation_WH[[1]][["coefficients"]][1,1]*real_wife$wife_age+
reservation_WH[[1]][["coefficients"]][2,1]*real_wife$wife_inc/10000+
reservation_WH[[1]][["coefficients"]][3,1]*real_wife$wife_ysh
real_wife$WH_std[real_wife$response==0] =
reservation_WH[[2]][["coefficients"]][1,1]*real_wife$wife_age+
reservation_WH[[2]][["coefficients"]][2,1]*real_wife$wife_inc/10000+
reservation_WH[[2]][["coefficients"]][3,1]*real_wife$wife_ysh
#Mid 和 Fid 一样的是一家子
real_husband=real_husband[order(real_husband$Mid),]
real_wife=real_wife[order(real_wife$Fid),]
cor(real_husband$husband_age, real_wife$wife_age, use = "complete.obs")
cor(real_husband$husband_inc, real_wife$wife_inc, use = "complete.obs")
cor(real_husband$husband_ysh, real_wife$wife_ysh, use = "complete.obs")
cor(married_sample_di$rage, married_sample_di$sage, use = "complete.obs")
cor(married_sample_di$rinc,married_sample_di$sinc,use = "complete.obs")
cor(married_sample_di$rysh, married_sample_di$sysh, use = "complete.obs")
#年龄差
mean(real_husband$husband_age - real_wife$wife_age, na.rm = TRUE)
sd(real_husband$husband_age - real_wife$wife_age, na.rm = TRUE)
#收入差
mean(real_husband$husband_inc - real_wife$wife_inc, na.rm = TRUE)
sd(real_husband$husband_inc - real_wife$wife_inc, na.rm = TRUE)
#教育差
mean(real_husband$husband_ysh - real_wife$wife_ysh, na.rm = TRUE)
sd(real_husband$husband_ysh - real_wife$wife_ysh, na.rm = TRUE)
#真实匹配相关系数
cor(real_husband$husband_age, real_wife$wife_age, use = "complete.obs")
cor(real_husband$husband_inc, real_wife$wife_inc, use = "complete.obs")
cor(real_husband$husband_ysh, real_wife$wife_ysh, use = "complete.obs")
#年龄差
mean(real_husband$husband_age - real_wife$wife_age, na.rm = TRUE)
sd(real_husband$husband_age - real_wife$wife_age, na.rm = TRUE)
#收入差
mean(real_husband$husband_inc - real_wife$wife_inc, na.rm = TRUE)
sd(real_husband$husband_inc - real_wife$wife_inc, na.rm = TRUE)
#教育差
mean(real_husband$husband_ysh - real_wife$wife_ysh, na.rm = TRUE)
sd(real_husband$husband_ysh - real_wife$wife_ysh, na.rm = TRUE)
