5.1. Integrating external data sources

One of the advantages of adopting this structured approach is that it enables the knowledge base we have created to include links to other external data sources. This enables us to support queries that seamlessly combine information from the impact case studies and these external data sources. For example, when a company appears in an impact case study, we create a link to its entry in DBpediaFootnote ¹ and to its Companies House records.Footnote ² This provides additional information, including the Standard Industrial Classification Code. An example is that we can answer the question: “Where are the companies included in impact case studies based?” (shown in Supplementary Appendix 7). The answer is shown in Table 7.

Table 7. Countries in which companies included in impact case studies are based, as measured by the count of unique company entities in “Unnderpinning research,” “Details of the impact,” and “Sources to corroborate the impact” found in Sections 2, 4, and 5 of the case studies

It would not be possible to answer this question if we could only use information held in the impact case studies. As can be seen from the country-level visualization presented in Figure 6, the research in UoA-11 has industry links to 39 countries, 38% of the companies were in the UK, followed by 24% in the United States, and 3.6% in the BRIC (Brazil, Russia, India, and China) nations.

Figure 6. The global reach of industrial impacts arising from research undertaken in UoA-11.

6.1. Assigning grades to individual case studies

Supervised machine learning algorithms require a training set of examples that include both the input features and the output (labels) that the system is trying to predict—in this case, the grade awarded to the case study. This presented a problem, as the grades for individual impact case studies were not made public. Instead, each case study was given a grade on a 9-point scale comprising integer and half-integer scores from 0 to 4. This was not made public. However, the grades awarded to all the Impact Case Studies in a Unit of Assessment were then combined, together with an impact template (described below) to give a single overall distribution that was published.

The distribution gives the percentage of case studies that were “unclassified,” 1*, 2*, 3*, and 4* (a case study or impact template given half-integer score had half of its grade assigned to both of the grades that it fell between). Each Unit of Assessment also submitted one impact template that described its approach to supporting and enabling impact from research conducted within the unit. The impact template was also graded on the same scale as the case studies, and contributed 20% to the overall distribution.

For example, if there were two case studies in the submission, each would contribute 40% to the submission, and the impact template provided the remaining 20%. If one case study scored 2.5, it would contribute 20% to both 2* and 3*. If another case study was graded as 3* and the impact template was graded as 1*, this would result in the overall published distribution: 1*: 20%, 2*: 20%, 3*: 60%, and 4*: 0%.

This scoring scheme makes it impossible to recreate the actual grade for every individual impact case study, which is needed for machine learning. To overcome this, the learning phase only used case studies from institutions where the variance in the grades awarded across all the case studies was low. This approach allowed us to use the weighted score (equation (1)) of the Unit of Assessment’s overall distribution as a reasonable estimate of each case study’s grade in a Unit of Assessment.

(1) $$ \mathrm{Weighted}\ \mathrm{score}=0.04.\left(4\ast \right)+0.03.\left(3\ast \right)+0.02.\left(2\ast \right)+0.01.\left(1\ast \right). $$

6.2. Setting the variance threshold

This approach raises the question of how to select the level of variance below which a unit of assessment’s case studies should be used by the machine learning algorithm. If the threshold is too low, then too few case studies will be eligible to train a robust machine learning model, whereas if it is too high, then the grades assigned to the individual case studies will not be accurate enough to train a reliable predictor.

In order to set a threshold, experiments were conducted across a range of 42 different variance thresholds. For each experiment, only the case studies from those institutions whose variance was below the threshold were used to train the machine learning model.

Figure 7 shows each predictor’s cross-validated accuracy (y-axis) using different numbers of training sets (x-axis) depending on the variance threshold. The labels on the x-axis show that as the variance threshold rises, the impact profiles from a greater number of institutions can be included in the machine learning. The graph shows that, after an initial sharp rise, there is a downward trend in the accuracy as the size of the training set increases—the accuracy ranges from a maximum of about 85% at 24 to a minimum of −140% at 72. A rapid decline in model performance starts at 47, where the predictor can still achieve acceptable accuracy. This was therefore the threshold we selected; these 47 impact profiles contain a total of 125 case studies.

Figure 7. The impact of using different variance thresholds to predict scores.

We then explored the distribution of grades of the impact case studies in the selected low variance profile. Almost all classification learning methods share an underlying assumption that the number of training samples in each different category is roughly equal. Imbalanced data may cause biases in the training process as it will be possible to achieve high accuracy by biasing predictions toward the major classes (Chawla et al., Reference Chawla, Japkowicz and Kotcz2004).

Figure 8 shows the imbalance across different grades for the impact case studies in the low variance profile (in orange).

Figure 8. Grade distribution of different parts of the impact case studies.

We used resampling methods to minimize the negative effects of imbalance. These methods modify the training dataset so that standard learning algorithms can be effectively trained on it. We, therefore, used ADASYN (He et al., Reference He, Bai, Garcia and Li2008), an adaptive synthetic sampling method for imbalanced data, to increase the size of the minority classes in our training set, choosing the sampling strategy based on the distribution given by the REF overview report (Higher Education Funding Council for England, 2015) (the blue line in Figure 8). We also applied “RepeatedStratifiedKfold” cross-validation during the oversampling, as (Santos et al., Reference Santos, Soares, Abreu, Araujo and Santos2018) suggest. The distribution of data after resampling is shown by the green line in the same figure.

Out of a total of 253 case studies submitted in Unit of Assessment 11, our approach allowed us to use 125 case studies for machine learning. The labels produced in this way are real numbers, which are only suitable for regression-based machine learning. To open the opportunity to explore classification-based machine learning algorithms we also created two other sets of integer labels:

• • 5-class labels—integers ranging from 0 to 4 (achieved by rounding the real number labels to the nearest integer)

• • 9-class labels—integers ranging from 0 to 8 (achieved by doubling the real number labels before rounding them to the nearest integer)

6.3. Feature selection

To accommodate the needs of machine learning, the key parts of the data are flattened from the RDF graph into a tabular format. The conversion from the graph to flattened data used Protégé (Musen and Protégé Team, Reference Musen2015) first in order to reduce human effort as much as possible. However, human input was applied to remove some very rare features that only exist in one or two specific case studies, as those features are too sparse for traditional machine learning tasks. The 42 flattened features used as inputs to the machine learning model fall into two categories: data from the whole of the unit of assessment (Table 8) and data from the impact case study itself (Table 9).

Table 8. Unit of assessment features

Table 9. Case study features

Having a large number of features increases the computational cost, and the risk of overfitting, while the number of observations required to produce a convincing prediction grows exponentially with the number of features. To prevent this, feature selection and extraction techniques are used to avoid this curse of dimensionality (McLachlan, Reference McLachlan2004). We employed only feature selection in this experiment since it preserves the original data characteristics when a subset of the features is selected. While feature extraction techniques such as Principal Component Analysis (PCA) are highly effective at reducing dimensionality, interpretability is necessarily reduced because the retrieved features are transferred to a new space. Interpretability was critical for this study since we wanted to determine which characteristics are critical for grade prediction.

Calculating all possible feature combinations prior to deleting unnecessary and redundant features is computationally impractical. As a result, we used the Deterministic Wrappers approach in this experiment, which follows the sequential forward selection (SFS) principle by greedily evaluating all alternatives, to match the 5-classes label based on an XGBoost classifier. The experiments were conducted individually for two categories, with the results displayed in Figure 9. Figure 9a illustrates the high accuracy achieved when five UoA features are employed; Figure 9b illustrates the initial increase in performance to 11 Case Study features before the accuracy begins to fall due to the number of dimensions.

Figure 9. Feature selection using the deterministic wrappers method.

To replicate the environment in which the REF panel made their decisions, we excluded those features (“number of degrees awarded,” “number of outputs,” and “GPA of outputs”) that were not available to the panel when they chose the grades.

Apart from the features described above, the symbol ^† is used in the glossary tables to denote other features not used in machine learning. The feature selection procedure indicated that these features only slightly improved the model’s performance, while increasing its dimensionality. The reason why some of these features do not contribute to the overall performance as much as might have been expected is that either some of them are highly correlated with each other, or they are Boolean features that have low variance distributions.

A clustered heat map (Figure 10) produced by the Seaborn package (Waskom, Reference Waskom2021) was used to visualize correlation analysis. As the features comprise both continuous and discrete variables, we used Spearman’s rank correlation coefficient to correlate the remaining features after feature selection and the estimated scores. Each intersecting cell shows the correlations between the variables on each axis. The strength of the color within its color gradient is proportional to the correlation coefficient. The dendrogram reveals the similarities between the variables. For example, the features “GovUse” and “Policy Influence,” the “SpinOff” and “Acquisition” were close in the case study context.

Figure 10. Heatmap visualization for selected features.

Rich information can be drawn from both positive and negative correlations between two variables in the heatmap. The main message is that “Research Income” is the variable most correlated with the score. The rest of the variables mostly follow the pattern: the higher the correlation with “Research Income,” the higher the correlation with the score.

6.4. Model performance metrics

The ability of the regressor and classifier to generalize from the test set and make accurate predictions of the grade are evaluated using the following commonly used evaluation criteria (Bishop, Reference Bishop2006) for supervised machine learning regression and classification models:

• • R-squared (R $ {}^2 $ ): The coefficient of determination measures a regression model with two sums of squares: one is the total sum of squares (denoted as SS $ {}_{\mathrm{TOT}} $ ) and the residual sum of squares (denoted as SS $ {}_{\mathrm{RES}} $ ). The definition is:

  (2) $$ {R}^2\hskip0.35em =\hskip0.35em 1-\frac{{\mathrm{SS}}_{\mathrm{RES}}}{{\mathrm{SS}}_{\mathrm{TOT}}}\hskip0.35em =\hskip0.35em 1-\frac{\sum \limits_i{\left({y}_i-{\hat{y}}_i\right)}^2}{\sum \limits_i{\left({y}_i-\overline{y}\right)}^2}. $$

• • Mean absolute error (MAE): used to show the average absolute error between the true value x $ {}_i $ and predicted value y $ {}_i $ with the same scale.

  The confusion matrix shown in Table 10 helps in understanding the classification metrics below.

  Table 10. The confusion matrix of classification

  

• • Precision is the probability of the sample actually being positive out of all the samples that were predicted to be positive.

• • Recall is the probability of a positive sample being predicted in a sample that is actually positive.

• • The F1-Score: In general, precision and recall are inversely proportional, we can either use the area under the Precision-Recall curve or F1-Score, which is the harmonic mean of Precision and Recall, to evaluate the classifiers from both perspectives. The formula is:

  (3) $$ \mathrm{F}1\hskip0.35em =\hskip0.35em \frac{2}{{\mathrm{Recall}}^{-1}+{\mathrm{Precision}}^{-1}}\hskip0.35em =\hskip0.35em 2\cdot \frac{\mathrm{Precision}\cdot \mathrm{Recall}}{\mathrm{Precision}+\mathrm{Recall}}\hskip0.35em =\hskip0.35em \frac{\mathrm{TP}}{\mathrm{TP}+\frac{1}{2}\left(\mathrm{FP}+\mathrm{FN}\right)}. $$

• • Hamming Loss: a loss function that indicates the fraction of the wrong class to the total number of classes, the optimal value is zero, a smaller value means better performance of the classifier (Tsoumakas and Katakis, Reference Tsoumakas and Katakis2007).

• • Area under ROC (AUC): defined as the area under the receiver operating curve (ROC) (Hand and Till, Reference Hand and Till2001). A ROC curve plots two parameters: sensitivity (a synonym for recall) and specificity.

The metrics above could be used directly in binary classification problems, whereas in multiclass classification problems, precision and recall should be averaged in the manner of Macro, Weight, and so forth (Sokolova and Lapalme, Reference Sokolova and Lapalme2009). Macro averages the metrics (Precision/Recall/F1-score) of the different categories, giving the same weight to all categories. It allows each category to be treated equally. The weight method gives different weights to different classes (the weights are determined according to the proportion of class distribution), and each class is multiplied by the weights and then summed. This method considers the imbalance of the classes and its values are more likely to be influenced by the majority class, which are grade 3, 3.5, and 4 in our case.

6.5. Predicting grades of impact profiles

We mainly used XGBoost (Chen and Guestrin, Reference Chen and Guestrin2016) and the built-in algorithms in scikit-learn (Pedregosa et al., Reference Pedregosa, Varoquaux, Gramfort, Michel, Thirion, Grisel, Blondel, Prettenhofer, Weiss, Dubourg, Vanderplas, Passos, Cournapeau, Brucher, Perrot and Duchesnay2011) for machine learning predictions.

XGBoost is a decision-tree-based ensemble machine learning algorithm that uses the gradient boosting framework. We chose it because it gives excellent results with small data sets, which we have in this situation. Additionally, we compared its performance to that of a Support Vector Machine (SVM) model.

The first experiment was conducted using the average impact grades of each institution in UoA-11 as the target label, with each university’s research income, and full-time equivalent (FTE) staff used as input features for training the model. This is to explore whether it is possible to use UoA features to make grade predictions over the impact profiles, which have a much smaller size compared to the number of impact case studies.

Comparison is made between two different subsets of the impact profiles. One is the full profile of all 89 institutions (the red legend in Figure 8), while the other is the 46 profiles that fall within the low variance threshold (the orange legend in Figure 8).

The experiments used 75% of the observations for training, and the remainder for testing. The results are shown in Table 11. In the XGBoost based regression model, the larger dataset only achieves 35% accuracy, while the smaller one achieved a higher accuracy (47%). This experiment shows the limits of using only the limited UoA features for impact profile grade prediction. It also shows the advantage of training on the low variance dataset, and so that is used for the rest of our experiments. To validate the performance of the case studies from the selected low variance impact profiles, we used the learning curves of three classifiers in Figure 11 to show the increase in macro-averaged F1-score performance of three feature sets—case study features, UoA features, and both of them when using different proportions of training data.

Table 11. Results for the full and low variance profiles using only UoA features

Figure 11. The learning curve of three different feature sets using the case studies within the low variance threshold.

6.6. Predicting the grade of a case study within the low variance profile

In the second experiment, two algorithms are applied to predict the scores on three different subsets of features using regression and classification methods. For the classification problem, the “5-classes” and “9-classes” labels defined above are used as the target for prediction.

For the classification experiment, we used 60% of the whole dataset, which is 75 out of 125 observations, for training and validation; the remaining 50 are for testing. We chose 60% over the commonly used 75% for the training set because there are only six samples for minor classes from 1 to 2. Increasing the test set to 40% allows the minor classes sample to increase to 18, making the minor classes evaluation less likely to be affected by having only a few samples. The learning curve in Figure 11 also shows the reasonable performance of three different feature sets when using 75 case studies.

Tables 12 and 13 show results from the nine experiments using XGBoost and SVM for regression and classification. We noted that the figures only show four classes (Figure 12); this is due to the lack of case studies with label 0 using our estimation method. Similarly, in the 9-classes classification, the figure only shows seven classes because of the lack of data with labels 0 and 0.5.

Table 12. XGBoost and SVM’s classification performance (Precision Recall and F1-Score) with 5 and 9 classes. The bold numbers indicate the best result for each feature set in weighted and macro average.

Table 13. XGBoost and SVM regression performance. The bold numbers indicate the best result for each feature set.

Figure 12. Classification report using XGBoost with three different features.

6.7. Model validation

We used classification reports to observe the model performance for each category through a confusion matrix, which comprises True Positives, False Positives, True Negatives, and False Negatives. Each classification report (Figure 12) shows the model’s performance using the metrics (Precision and Recall and F1-Score) corresponding to each category. Two other metrics (AUC (Area Under Curve) and Hamming loss) for XGBoost classifiers are provided in Table 14.

Table 14. XGBoost AUC and Hamming loss with 5 and 9 classes. The bold numbers indicate the best result for each feature set.

From those results listed above, we can see the XGBoost model outperforms the SVM in our dataset, so we will mainly focus on interpreting the results from XGBoost models. The UoA features achieve satisfactory performance in Grade 3–4, but much lower performance in Grade 1.5–2.5. The case study features achieve relatively lower performance compared to the UoA features in Grade 3–4, but have a more stable results across all seven grades. This may explain why the UoA features outperform case study features in the 5-classes classification but have similar performance with case study features in the 9-classes classification. When we use the combination of the two features, the performance in Grade 3–4 remains as high as the UoA features and outperforms UoA and case study features separately. Different metrics using AUC and Hamming loss in Table 14 also tell the same story. We also used the validation curves to show the generalization behavior of the classifiers using UoA and all features in the 9-classes classification (Figure 13). Regardless of the two small gaps between the two validation curves, which indicate the need for a larger data set, neither the UoA classifier nor the all features classifier has under-fitted or over-fitted.

Figure 13. Validation curves of classifiers using UoA and all features in 9-classes classification.

6.8. Result interpretation

Several methods are provided in the XGBoost package to interpret the model that has been created. Figure 14 provides a visualization of part of one of the trees in the model that split the data using the conditions. However, it is not practical to visualize the full trees given their size.

Figure 14. Tree visualization of the decision tree for grading case studies.

Therefore, a different method was used to measure the feature importance for the model shown in Figure 12c, which uses both UoA and case study features as input and 5-classes as output. Figure 15 shows feature-importance ranked by weight, coverage, and total gain:

• • Weight is the frequency of a specific feature used for splitting all trees across all weak learners. This shows that the “Reference Cite” feature has been used the most in splitting the trees.

• • Coverage is the number of observations that are split, weighted by the frequencies of the same feature used in splitting across all the trees. This indicates that the feature “Patent” has split the most observations.

• • Total gain is the total gain in accuracy brought by the specific feature when splitting the trees. This shows that the feature “Research Income” contributed the most to the model prediction accuracy.

  

  Figure 15. Feature importance (Weight, coverage, and total gain).

It should be noted that binary features like “Improved Environmental Sustainability” can only be used once in each tree. As a result, binary features are scored much lower in metrics where Weight (frequency) is used.

For a more in-depth analysis of the results, we use SHAP (SHapley Additive exPlanations) (Lundberg and Lee, Reference Lundberg and Lee2017) to provide a more consistent and accurate way of measuring features. SHAP is an innovative approach to improving the predictions from a complex model and exploring relationships between the features and the individual observations.

Figure 16 shows a force plot produced by SHAP values. The blue side indicates a negative contribution to the prediction, while the red side indicates a positive contribution. As shown on the right-hand side, the features “Research Income,” “Collaboration,” “Duration,” “Open Source,” and “FTE” have a positive contribution, while the “Total Funding in GBP” and “Academic Use” have a negative contribution. This may be due to the fact that a significant proportion of the impact case studies do not mention financial information or only mention very vague figures. Of the 125 observations, 24 out of 31 case studies had blank funding data and scored between 3 and 4. Interestingly, the remaining seven cases only scored between 1 and 1.5. This lack of information may have led to features such as “Total Funding in GBP” being ranked low in terms of feature importance.

Figure 16. SHAP force plot.

To get an overview of which features are most important for a model we can plot the SHAP values of every feature for every sample. The summary plot in Figure 17 sorts features by the sum of SHAP value importance over the case studies and uses SHAP values to show the distribution of the impacts each feature has on the model output. The color represents the score value (red high, blue low). This reveals that high “Research Income,” “FTE,” “Duration,” “IndustrialUse,” “SpinOff,” “Standard,” “PublicUse,” “Acquistion,” and “AcademicUse” increases the predicted score; higher “KTP” lowers the predicted score.

Figure 17. SHAP summary plot.