R Notebook

library(dplyr)
library(matchingMarkets)
library(ggplot2)

result_vector <- c(NULL)
for(times in 1:1000){
  student_random_n <- matrix(rep(0,400),nrow = 4)
  for(i in 1:100){
    student_random_n[,i] <- sample(1:10,4)
  }
  
  school_random <- matrix(rep(0,1000),nrow = 100)
  for(i in 1:10){
    school_random[,i] <- sample(1:100,100)
  }
  #match matrix
  match2 <- matrix(rep(0,1000),nrow = 100)
  for(i in 1:100){
    for(j in 1:4){
      k <- student_random_n[j,i]
      match2[i,k] <- 1
    }
  }
  
  for(i in 1:10){
    n=0
    x=1
    rank <- school_random[,i]
    for(j in rank){
      x=x+1
      if(match2[j,i]==1){
        n=n+1
      }
      if(n>=20){
        rank2 <- rank[x:100]
        for(k in rank2){
          match2[k,i] <- 0
        }
        break
      }
    }
  }
  
  #renew preference
  student_random_n2 <- student_random_n
  school_random2 <- school_random
  
  for(i in 1:100){
    for(j in 1:4){
      k <- student_random_n[j,i]
      if(match2[i,k]==0){
        student_random_n2[j,i] <- 0
      }
    }
  }
  for(i in 1:10){
    for(j in 1:100){
      k <- school_random[j,i]
      if(match2[k,i]==0){
        school_random2[j,i] <- 0
      }
    }
  }
  for(i in 1:100){
    for(j in 1:4){
      if(student_random_n2[j,i]==0){
        student_random_n2[,i] <- c(student_random_n2[,i][-j],0)
      }
    }
  }
  
  for(i in 1:10){
    k <- c(NULL)
    for(j in 1:100){
      if(school_random2[j,i]==0){
        k <- c(k,j)
      }
    }
    a <- school_random2[,i][-k]
    b <- rep(0,length(k))
    school_random2[,i] <- c(a,b)
  }
  
  #DA Algorithm
  res_r_r <- hri(s.prefs = student_random_n2,c.prefs = school_random2,nSlots = rep(10,10))
  res2_r_r <- res_r_r$matchings
  res2_r_r <- as.data.frame(filter(res2_r_r,sOptimal==1))
  result <- length(res2_r_r$matching)/length(unique(res2_r_r$matching))
  result_vector[times] <- result
  if(times%%100==0){
    print(times)
  }
}

[1] 100
[1] 200
[1] 300
[1] 400
[1] 500
[1] 600
[1] 700
[1] 800
[1] 900
[1] 1000

resvdf <- data.frame(result=result_vector)
resvdf$times <- c(1:1000)
ggplot(resvdf,aes(result))+
  geom_bar()+
  scale_x_continuous(breaks = c(85:100))+
  geom_vline(aes(xintercept=mean(result)),color="Black", size=1)+
  theme_light()