Work Visa: Ensemble Methods¶

The number of applications for US work visas is growing. To assist with the review process, we will develop an ML solution to filter out many candidates. This filtering will be trained on a data set of previous applications and the resulting case statuses.

This data has features describing both the applicant and their sponsoring employer:

• case_id: unique identifier for each application

• continent: employee continent of origin

• education_of_employee: level of education

• has_job_experience: binary flag

• requires_job_training: binary flag

• no_of_employees: size of sponsoring employer's company

• yr_of_estab: sponsoring company's year of establishment

• region_of_employment: applicant's intended region of employment in the US

• prevailing_wage: regional average wage for a given domain of labor; useful to keep job market competative without underpaying foreign workers

• unit_of_wage: frequency of payment

• full_time_position: binary flag; Y: full time

• case_status: binary flag

Importing necessary libraries and data¶

Data Overview¶

There are twelve features with 25480 records each. While we record no NaN values, there may very well be missing or erroneous data. Most of the columns are object type, which we can convert to categorical during our pre-processing later.

Exploratory Data Analysis (EDA)¶

• We see that the minimum number of employees is -26, an impossibility! This indicates that we already have at least one apparent error in our data.

• The median number of employees is around 2100, so a midsize company. The largest in our records has over 600,000 employees.

• The oldest company we have records for was founded in 1800, but the mean year is around 1979.

• The median local wage is $70,308.21, and the mean is several thousand higher, indicating some skew toward higher costs of living. There is at least one shockingly low value, as the minimum is just $2.14. Perhaps this is hourly wage, though even then $2/hr is quite low. This hints that we must diligently explore outliers for potential errors.

• With 25480 unique case IDs, we can be fairly confident that there are no duplicates in our records.

• The data covers 6 continents of origin, with Asia being the most frequent.

• The plurality of applicants have a Bachelor's level of education. We will explore the other levels on record shortly.

• The majority of applicants do have some job experience, and around 88% of jobs do not require training.

• Nearly all wages are recorded as yearly salaries, though the data dictionary indicates that hourly, weekly, and monthly are options too.

• The only continent of origin not represented is (predictably) Antarctica. After Asia, the most common continents are Europe and North America (likely Canada and Mexico).

• Other than a Bachelor's degree, many applicants have a Master's. Some also have a high school diploma or a doctorate, but these are less common.

• There are five options for region of employment: Northeast, South, West, Midwest, and Island (in decreasing level of frequency). Northeast and South have nearly the same number of records.

• The dependent variable in this study is case_status, with outcomes Certified or Denied. About two-thirds of cases are Certified.

We find that Asia is by far the most frequent continent of origin, while Oceania is the least.

Bachelor's and Master's education are almost equally common, with both around 10000 records each.

The majority of applicants have some job experience, but over 10000 (about 40%) have none.

Unlike the previous feature, which had fairly balanced classes, most applicants do NOT require job training.

With such skewed data, the histogram lends little insight. Boxplots give a better idea of data concentration and distribution. Most of our records are concentrated below around 7000, but there are numerous extreme values on the high end.

We will not treat these extreme values though. It is entirely reasonable that many records reflect applicants to large companies.

We find that 75% of companies in our records were founded in 1976 or later. Like the previous feature, we have many extreme values, this time on the lower end. As before, all these values are entirely sensible, so we will not alter them.

However, we will convert this feature to 'years since founded', which will flip the distribution from left-skewed to right-skewed but otherwise leave the data unchanged.

Northeast, South, and West are all common regions. Island is decidedly uncommon, with fewer than 500 records.

Aside from the spike around $0, prevailing wage follows a right-skewed normal distribution. It may well be that the spike is due to other wage units, such as hourly.

Indeed, the distribution looks far less erratic when just focused on yearly prevailing wage.

Most of our records consider cases with yearly wages. Weekly and montly are vastly less common, as was shown earlier in the value counts table.

Most positions are full time.

Only one-thirds of cases are Denied.

Let's dig deeper into prevailing wage.

I am skeptical about the records tagged with 'Week' for unit_of_wage. According to our records, 75% of jobs with a weekly wage pay at least $51,000 per WEEK! As we have limited visibility into the source of these data, we will preserve these records, but ideally I would like to do further research to justify these entries.

The boxplot above is further evidence that weekly—and even monthly—data looks erroneous. Many of the weekly wages would result in a yearly earning well over a million USD! Perhaps this is personal bias, but it seems far fetched to have this many exceptionally high-paying jobs in this data set.

Without any evidence to the contrary, however, we will leave these records intact.

For the record, some hourly data points seem unreasonable too: The maximum wage is about $1000 per hour, or around $2 million per year. Again, there is no evidence that this data is necessarily erroneous; it simply stands out.

The middle 50% of the data for prevailing wage is higher in the midwest and island regions. The northeast has the lowest first quartile. Every region has many extreme values on the high end.

Looking toward our dependent variable, it appears that case_status is influenced by education. The trouble here is that it is difficult to compare the regions, as the scale is different for each. Is the percentage Certified for Doctorate any greater or less than that for, say, Bachelor's?

Let's instead examine the percentage Certified and Denied.

Indeed, here we can compare case_status on like scales. We find that applicants with a Doctorate are most often certified (over 87% of the time), while employees with only a High School education are certified around a third of the time. We find a clear trend: Visa status is directly connected with level of education. Higher levels of education lead to a greater percentage of visas certifed.

To plot the rest of our features as percentages, we make the above code into a function.

• Europe, Africa, and Asia see the greatest percentage of certified visas.
• South America sees the least.
• Every region has over 50% certification rate.

Having job experince certainly helps, as the proportion of Certified visas is higher for applicants with job experience. Let's test whether this difference is significant, assuming a level of siginificance of 5%. Please note that our populations are binomially distributed, independent, and simply randomly sampled; they furthermore satisfy the sample size inequalities. The assumptions for a two proportion z-test are therefore met. Our null hypothesis is that our sample proportions come from populations with the same ratios of Certification to Denial; the alternative is that the ratios are truly different.

With an astoundingly low p-value, we can confidently conclude that these sample proportions reflect a real-world difference for visa applications: We find that a greater proportion of applications are certified when the applicant has previous job experience.

Unlike the previous feature, the requirement of job training has very close proportions: 66.6% Certified (N) versus 67.9% (Y). There is apparently no difference in case_status based on the job training requirement. To confirm, we can run another two proportion z-test. (As before, the assumptions for the test are met.) We set the level of siginificance at 0.05. Our null hypothesis is that our sample proportions come from populations with the same ratios of Certification to Denial; the alternative is that the ratios are truly different.

Our p-value exceeds 0.05, so we fail to reject the null hypothesis. We must conclude that there is no significant difference between the two proportions. Put another way, requires_job_training does not impact case_status on its own.

There seems to be little difference in visa certification rate among island, western, and northeastern jobs. The south and midwest look to have appreciably higher rates.

Applicants to hourly jobs see far less than 50% certification. This is well below the rest, with yearly salaried jobs being most likely to lead to visa certification.

Similar to requires_job_training, it appears here that there is not much difference in the resulting certification based on whether the position is full time. We find that 68.5% are Certified for part time work and 66.6% are Certified for full time work. Our assumptions being once more satisfied, we can run a two proportion z-test. We set the level of siginificance at 0.05. Our null hypothesis is that our sample proportions come from populations with the same ratios of Certification to Denial; the alternative is that the ratios are truly different.

With a p-value under 0.05, we can reject the null hypothesis. In this case, there is a statistically significant difference between these two proportions. That means that an applicant is more likely to be certified for part time work than full time work.

The mean prevailing wage for certified applications is higher than that for denied applications by about $8500. The middle 50% for certified is concentrated higher too. Interestingly, the whiskers extend further for denied applications. Taken with its wider IQR, this shows that denied applications have greater variance in prevailing wage.

This is also true across various regions.

Denied applications in the northeast seem to have some of the lowest prevailing wages. The prevailing wages for the island region, on the other hand, appear generally higher than the rest, with certified cases having the greatest middle 50% among the regions. This is likely due to the higher cost of living on islands, as most goods must be imported.

Data Preprocessing¶

Feature Engineering¶

Rather than recording the year businesses were established, we will instead consider the number of years in business.

Outlier detection and treatment¶

There are 33 records where the application lists a negative number of employees. These clearly constitute an error, and with no clear method to impute these values, our best recourse is to drop the rows.

Data Prep¶

We drop case_id, as the unique identifiers will not lend any predictive power to our model.

After making all object variables into category type, we encode them numerically.

We use one hot encoding for the categorical variables with no inherent order.

Now that we've split our data, it is ready for modeling.

Second EDA¶

We only removed 33 records out of 25480. That's only 0.1% of the data. Accordingly, there is not enough change in our data to produce appreciably different visualizations. We only include plots for features that were directly altered/engineered.

• We find that most businesses in our records are less than 50 years old.
• There are a steady supply of businesses greater than 75 years old, about the same amount for each year. It only starts to trail off after around 175 years.

With so few records removed, there is little difference in the plots for number of employees. Put another way, the comments on this feature from the first EDA are still applicable. The middle 50% of observations is between about 1000 and 3500. There are also tiny companies, with not more than several employees, and massive companies, with over half a million employees, in our data set.

Building bagging and boosting models¶

When scoring our models, we must find a good balance between false positives and false negatives.

• A false positive is a certified visa that should have been denied. We want to reduce false positives because these workers may be ineffective and not benefit the company.
• On the other hand, a false negative constitutes a denied visa that should have been certified. This deprives the workforce of an asset, someone who could benefit US companies.

With this in mind, we will score our models on F1, as it balances recall and precision, reducing both false positives and false negatives.

Functions¶

This function gets the scores of a model on training data and testing data.

This function adds accuracy and F1 scores to a DataFrame that we can use to compare models. It also runs a two proportion z-test to determine whether the training accuracy and testing accuracy are different enough for the model to be considered overfit. We assume a level of significance of 5% in the interpretation of the p-value. All conditions are necessarily met to justify the test.

The function prints a heatmap with the confusion matrix for a given model.

Bagging¶

Baseline: Decision Tree¶

As a baseline comparison, we train a single decision tree.

Predictably, this model is comically overfit.

We add this model to a DataFrame that we will use to compare the different models.

Bagging Classifier¶

The default estimator on a bagging classifier is a decision tree classifier, so this model consists of many parallel decision trees trained on samples taken with replacement (bootstrap samples).

This model is also overfit. The test performance is nonetheless slightly better than the single decision tree.

Random Forest Classifier¶

As with the others, the random forest classifier overfits on the training data. Test performance is better than before, though this model is not generalized.

Boosting¶

AdaBoost¶

The basic AdaBoost Classifier has good performance and does not suffer from overfitting. Accuracy on training and testing is around 73% and F1 is about 82%. This is a much more promising baseline.

Gradient Boosting¶

Gradient Boosting also does not appear overfit, with the added bonus of slightly better performance over AdaBoost.

• Train and test accuracy is around 75%.
• Train and test F1 is around 82%.

XGBoost¶

XGBoost offers near identical performance to the gradient boosting above.

Will tuning the hyperparameters improve the model performance?¶

Bagging¶

Tuned Decision Tree¶

• We will test values of tree depth from 3 to 9.
• We will consider limiting the maximum number of features to consider when calculating the best split.
• We will try balancing the class weights, though I don't expect this will improve performance.

Grid Search Findings

• Our response variable is not so unbalanced that a rebalancing of class weights is required.
• A maximum tree depth of 7 seems to be ideal, with 20 as the upper bound for leaf nodes.
• The model requires at least 5 samples for each leaf. Along with limiting the depth and maximum number of leaves, this reduces overfitting.
• We get best performance when considering all features to determine the best split.

We refit the best estimator on the full training set. (In the grid search, we employ cross validation, so there is a holdout set used to test each model. The upside is a more reliable search, but the cost is that we have yet to train the best estimator on the full training set.)

With around 75% accuracy, the tuned decision tree has similar performance to untuned boosting models. F1 score is around 82%. Importantly, tuning this tree has eliminated overfitting.

Tuned Bagging Classifier¶

• We will test different depths for the base estimator (Decision Tree Classifier).
• We will consider between 10 and 50 estimators.
• Max samples and max features will also be adjusted, testing different fractions to include.

• Trees with a depth of 5 score best on F1.
• Our model performs better with fewer estimators.
• Taking less than 100% of the samples yields a higher-scoring model.

Note that the tuned bagging classifier has great recall. This model does a great job at reducing false negatives, meaning the US job market is not losing out on promising talent. On the other hand, precision is lower, which leaves us with an F1 score comparable to other generalized models.

As an experiment, we try another grid search with a second option for the base estimator: Logistic Regression.

After fitting the grid search, we find that the decision tree classifier was still the strongest base estimator for our metric. However, this time more estimators were required.

It appears our model now overfits, as we have not limited the tree depth. We will not consider this model.

Tuned Random Forest¶

• We consider a number of estimators between 150 and 300. We also consider a huge model: 1000 trees.
• Various depths will be studied.
• As with the previous bagging model, we will test various values for the maximum number of features and samples.

• Limiting the number of features at each split by the square root of the total number of features yielded better results than by using log base 2.
• A maximum tree depth of 5 worked best.
• We required 250 estimators in the highest scoring model.

Similar to the tuned bagging classifier above, we're getting over 90% recall with this model. Unfortunately, accuracy is still below 75% and F1 is no higher than before.

Boosting¶

Tuned AdaBoost¶

Next we'll train an AdaBoost classifier, using a decision tree as the base estimator. We will consider between 10 and 70 estimators within the ensemble, various depths of the trees, and several different learning rates.

With a depth of 3 and a slower learning rate, our best model required 50 estimators.

We get similar F1 performance to the tuned random forest, but better accuracy. Recall is lower than some other models.

Tuned Gradient Boosting¶

We will next tune the gradient booster with several different initialization conditions. We will also check the number of estimators and the maximum number of features when looking for the best split.

As with AdaBoost, our accuracy is over 75%. The F1 score does not seem to be able to improve above 82%.

Tuned XGBoost¶

We will try many parameters for XGBoost.

• The parameter eta is the learning rate for the model.
• The parameter gamma concern loss reduction when partitioning a leaf node on a tree.
• Several sampling parameters control how the model samples the training data during fitting.
• The last parameter (scale_pos_weight) can help with biased classes in the response variable.

• A lower learning rate ended up being better for F1.
• The classifier did not need to adjust class weights for better performance.
• While colsample_by_tree is less than 100%, both colsample_by_level and subsample are maxed out at 1.

After training the best estimator on the full training set, we find that performance is comparable—in fact, nearly identical—to the other tuned boosting models.