This third and final part of the project follows the EDA part, and intends to dig deeper into the findings made in that EDA. In the first section, we build our statistical model, which includes defining the model, assessing its goodness of fit with various methods, and eventually refining the model; and in the second section, the final model is validated in order to have a better picture of its predictive power on new data.

This data set has the most relevant columns of merged_df from the EDA part, and 780 rows since the rows corresponding to athletes from the countries with the smallest numbers of observations were excluded (cf EDA).

Section 1: Building the model

Some techniques such as random forests cannot handle character information, so character variables need to be converted to type factor. For our model, the explanatory variables we are particularly interested in are meetYear, event, gender, ageAtMeet and country. However, country is a character variable that has a lot of possible values, which will yield many predictors and make our model too complex:

There are 25 unique country names in the country column, which if converted to a factor variable as it is, will give 25 levels. With the four other explanatory variables, our model will have well over 30 predictors, which sounds like a lot given our sample size of only 780 observations. This might lead to overfitting, or just simply make the results of our statistical analysis difficult to interpret. One way we could try to address this problem is replacing country in our model with another factor variable with a much smaller number of levels. In the previous EDA part, we interpreted the number of observations by country as an indicator of the strength of that country’s swimming tradition. Let’s have a look at the number of observations per country:

As a reminder, each observation (row) originally contained the best individual event time and the best relay split time for the same event in a certain meet. The more observations a country has, the more athletes from that country swam both individual and relay events which is a good indicator of the strength of a national swimming tradition. Given the number of countries and of observations, we could think of giving that new variable 3 levels: limited swimming tradition, moderate swimming tradition and strong swimming tradition. When it comes to allocating countries to these three categories, we are going to take into account the position of each country compared to the median and quartiles of the distribution, as well as the actual number of observations:

Countries with 4 to 13 observations will belong to the limited tradition category, those with 22 to 42 observations to the moderate tradition category, and those with 52 and more observations to the strong tradition category.

Now, what is left is to convert the other character variables to factor variables as well. But before we do so, there is another variable we need to transform if we want it to have any significant predictive power. Indeed, if we try to apply our model to the athlete prospects for the 2024 Paris Olympics or next year to the swimmers who will take part in the 2025 Worlds in Singapore, the majorMeet values for each future predicted observation will be new values (say, paris_2024 and singapore_2025), so it does not really make sense to include majorMeet as it is right now as a factor variable in our model. However, it certainly makes more sense to decompose majorMeet into 2 variables accounting for the type of major meet, which we will call meetType, and time, which is meetYear. In the EDA part, we identified what looked like a downward trend for faster relay TIW over the years at World Championships, while for Olympic Games we lacked data to identify a real trend, but we still might be able to capture potential trends linked to the type of major meet by introducing meetType. This is what EDA is for!

And finally, we convert event and gender to factor as well:

We can now build our model, which consists of fasterRelay, our dependent variable, and the 6 explanatory variables below:

Our model has 11 predictors (constant + 1 + 4 + 1 + 1 + 2 + 1), which is much better than 36 had we not substituted country for swimTradition.

Our dependent variable being binary, training it with logistic regression and classification trees (random forests, gradient boosting) is more adequate than linear regression. The first step in our analysis will be training the model with logistic regression, and evaluate how well it fits the data with measures of goodness of fit and likelihood ratio tests, and study the deviance residuals to check for heteroskedasticity. Once we have a good grasp of the predictive capabilities of our model, we will then train and test it with different methods (logit, simple classification trees, random forests, and gradient boosting) by performing K-fold cross-validation, and evaluate the outcome by studying the area under the resulting ROC curves for each classification method.

Evaluating the goodness of fit with logistic regression

To begin we first train our model with the glm function:

Then we define the observed values and make predictions:

The results for logistic regression are obtained by running the command below which yields Table 1 below:

stargazer(rsltLog, summary=FALSE, intercept.bottom=FALSE, align=TRUE, no.space=TRUE, type="html")

	Dependent variable:
	fasterRelay
Constant	93.033*
	(54.348)
meetYear	-0.047*
	(0.027)
event100m Breaststroke	0.066
	(0.255)
event100m Butterfly	0.228
	(0.256)
event100m Freestyle	1.521***
	(0.240)
event200m Freestyle	0.893***
	(0.244)
gender2	-0.255
	(0.156)
ageAtMeet	0.062***
	(0.023)
swimTradition2	0.723**
	(0.290)
swimTradition3	0.611**
	(0.260)
meetTypeworlds	-0.423**
	(0.169)
Observations	780
Log Likelihood	-491.756
Akaike Inf. Crit.	1,005.512
Note:	*p<0.1; **p<0.05; ***p<0.01

Table 1

Looking at the results in Table 1, we can see that several coefficients have asterisks next to them, indicating that their impact on our dependent variable is statistically significant according to the levels specified in the footer note. These results are quite pleasantly surprising since most of them corroborate what we observed in the EDA part.

This first logistic regression analysis found some evidence of meetYear having a slight negative impact on faster relay TIW ($\beta_{meetYear}=-0.047$, $p<0.1$), although meetType seems to explain most of the variation in faster relay TIW related to majorMeet, as the impact of the World Championships meet type over Olympics on relay TIW is stronger and more significant ($\beta_{meetTypeworlds}=-0.423$, $p<0.05$). This would mean that the year in which a major meet takes place, or the passing of time itself does not explain the variation in faster relay TIW nearly as much as the type of this major meet, i.e. Worlds or Olympics. In other words, if a swimmer who retired in 2015 said “In my day, swimmers used to care much more about the national team’s success than they do nowadays!”, these results right here could prove them wrong with high (but not absolute) certainty.

Moving on, we also notice that specialising in the 100m Breaststroke and the 100m Butterfly does not significantly influence faster relay TIW compared to the 100m Backstroke ($\beta_{event100m Breaststroke}=0.066$, $p>0.1$; $\beta_{event100m Butterfly}=0.228$, $p>0.1$). However, as we observed in the EDA part, being a 100m Freestyle or 200m Freestyle swimmer has a much stronger positive and significant impact on relay TIW than being a 100m Backstroke swimmer ($\beta_{event100m Freestyle}=1.521$, $p<0.01$; $\beta_{event200m Freestyle}=0.893$, $p<0.01$).

As we suspected during the EDA, being a female top-level swimmer does not significantly influence the outcome of relay TIW compared to male top-level swimmers ($\beta_{gender2}=-0.255$, $p>0.1$).

We notice a positive and highly significant impact of the age of a swimmer at the time of the meet on faster relay TIW ($\beta_{ageAtMeet}=0.062$, $p<0.01$), which means that as we discovered during the EDA, older swimmers are more likely than younger ones to swim faster TIW for their team than in solo, although this factor is not very strong, or doesn’t greatly impact faster relay TIW.

Finally, it turns out that being from a country with a moderate or strong swimming tradition has a rather strong, positive and significant impact on faster relay TIW as opposed to a country with a more limited swimming tradition ($\beta_{swimTradition2}=0.723$, $p<0.05$; $\beta_{swimTradition3}=0.611$, $p<0.05$).

Measuring the goodness of fit of mdlMain

We can first make a quick likelihood ratio test to find out if our fitted model really helps predict the outcome of fasterRelay better than a null model such as $fasterRelay \ \ \sim \ \ 1$:

The deviance of our fitted model is smaller than the null deviance, which indicates that the former does a better job than the latter at explaining the variation in the data, but to find out if that difference is statistically significant we build an LRT statistic (which follows a $\chi^{2}$ distribution) with those two deviances and look at the left-tail probability which is close to zero (p-value: 0). In other words, our null hypothesis being $H_{0}:$ “The fitted model mdlMain does not predict the outcome of fasterRelay better than the null model $fasterRelay \ \ \sim \ \ 1$.”, our p-value close to 0 tells us that we can reject it.

Unlike with linear regression, logistic regression does not have a straightforward and unique goodness of fit measure. In this study we will use McFadden’s pseudo $R^{2}$ which measures the improvement of our fitted model mdlMain over the null model. As a reminder, McFadden’s pseudo $R^{2}$ is defined as: $$R^{2}_{McF} = 1 – \frac{loglik(fitted)}{loglik(null)}$$ , so we first need to calculate these log-likelihoods.

The adjusted version of McFadden’s pseudo $R^{2}$ penalises extra predictor variables. For the sake of comparison, the count of correctly classified positives Count.R2 was computed, the adjusted version comparing the number of correct classifications to the frequency of the dominant outcome. We see that our McFadden’s pseudo $R^{2}$ is slightly below 0.1, which indicates that our current model does not do an excellent job at fitting the data. This pseudo $R^{2}$ is not interpreted exactly like the linear regression $R^{2}$ in the sense that we rarely see pseudo $R^{2}$ values of over 0.5, according to McFadden himself a pseudo $R^{2}$ of between 0.2 and 0.4 is already indicative of excellent goodness of fit (Behavioural Travel Modelling, 1979). Ours being well below this range tells us that there is room for improvement, which is corroborated by a Count $R^{2}$ of only 0.651 (only 65% of positives or fasterRelay = 1 were correctly classified).

All the GVIFs (variance inflation factors) obtained hover above 1 which at first sight does not raise concern for multicollinearity. This is not surprising, as we did not suspect particularly strong association between any of our explanatory variables. However, the GVIF cannot detect more subtle differences between performances of the levels of a factor variable.

In order to confirm the non-significance of gender and meetYear which both have coefficients with a p-value above 0.05, we perform a $\chi^{2}$ test on each of them to check for the significance of their individual effect on the predictive power of our model:

Unsurprisingly, both tests yielded p-values well above 0.05 which is our significance level, therefore these two variables are candidates for exclusion from our model.

Next, using the emmeans() (estimated marginal means) method we will be able to discern how the 5 levels of event contribute to the prediction of fasterRelay, and whether the differences in predictive power between them are statistically significant with the aim to potentially group levels with similar predictive power and reduce the complexity of our model. If a red (left and/or right) arrow of a category overlaps with an arrow of another category, it means that these two categories are actually possibly one.

        Estimated marginal probabilities by event

What the pairs() function does is basically compare the estimated marginal probabilities (event.emp) of each event against the others, in pairs. There are five events, which gives 10 pairs. The p-values were already adjusted for our 10 pairwise comparisons using the Tukey method, so we can directly make a decision on which levels to group depending on the p-value of their pair. We decide to set the level of significance at $\alpha = 0.05$. Our null hypothesis here would be $H_{0}:$ “Event A and Event B equally contribute to the prediction of fasterRelay“, where Event A and B can be any two different events from the five factors of event. If the p-value for a pair is greater than our $\alpha$, we can conclude that there is not enough evidence that Event A’s and Event B’s contributions to the prediction of fasterRelay are significantly different. Conversely, if the p-value of a pair is below $\alpha$, we can reject the null hypothesis for that pair. Given an alpha of 0.05, we can reject the null hypothesis for all pairs except three: 100m Backstroke - 100m Breaststroke ($p=0.999$), 100m Backstroke - 100m Butterfly ($p=0.901$), and 100m Breaststroke - 100m Butterfly ($p=0.972$). These three events were all compared against eacher other, and all the p-values for these comparisons were very high. Therefore, along with the results previously obtained in Table 1, we can conclude with high certainty that the 100m Backstroke, 100m Breaststroke and 100m Butterfly levels of event all equally contribute to the prediction of fasterRelay and can therefore be grouped into one unique factor.

In order to stay consistent in the assessment of statistical significance throughout our analysis, we must evaluate all variables with the same $\alpha$ which we set to 0.05, so to conclude our ANOVA chi-square tests and pairwise comparison test, following that level of significance we decide to exclude gender and meetYear from our model, and to group the 100m Breaststroke and 100m Butterfly levels of event with the reference level 100m Backstroke. Our new model will therefore have 7 predictors, down from 11.

Building and assessing our new model mdlMain.R

We are first going to define our new restricted model (hence .R) following the results of the previous section:

The first three levels are grouped into a single 100m Breast_Back_Fly level:

The previous levels were successfully replaced with the new level. Let us now build our new restricted model:

We also define a second model, mdlMain.Rf, identical to the first one except for fasterRelay being converted to factor which will be needed for the random forests method in the next section.

We once again perform logistic regression on our new model:

	Dependent variable:
	fasterRelay
Constant	-2.523***
	(0.617)
event100m Freestyle	1.416***
	(0.190)
event200m Freestyle	0.799***
	(0.195)
ageAtMeet	0.066***
	(0.022)
swimTradition2	0.689**
	(0.287)
swimTradition3	0.579**
	(0.258)
meetTypeworlds	-0.443***
	(0.167)
Observations	780
Log Likelihood	-495.131
Akaike Inf. Crit.	1,004.263
Note:	*p<0.1; **p<0.05; ***p<0.01

Table 2

The logit coefficients for our new model can be seen in Table 2, and this time every coefficient is statistically significant at least at a level of 0.05. The coeffficients of our explanatory variables slightly changed, but the interpretation of their influence on fasterRelay stays the same.

Just like we did for our previous model, we are now going to assess the goodness of fit of mdlMain.R:

A quick likelihood ratio test confirms that our new model still brings significant improvement in the prediction of fasterRelay, however the deviance of mdlMain.R looks pretty similar to that of mdlMain. Let us now have a look at the pseudo $R^{2}$:

The pseudo $R^{2}$ ($R^{2}_{McF} = 0.074$) and the Count $R^{2}$ ($R^{2}_{Count} = 0.627$) are both almost the same as those of mdlMain, approximately 5% smaller. This means two important things: the restricted mdlMain.R model does not fit our data better than mdlMain, but these goodness of fit results also confirm that the predictors we excluded from our new model indeed did not play a significant role in predicting fasterRelay.

Before we move onto model validation, the last step in this section is to perform a quick residual analysis, which will check for heteroskedasticity in our model.

Deviance residuals analysis

Another way to look at the goodness of fit of our model is to analyse the deviance residuals of our model, which are basically the components of the (binomial) log-likelihood function but slightly modified to follow a $\chi^2$ distribution (if we sum the squares of the deviance residuals and then divide by -2, we get the log-likelihood function). Let us have a look at both the Pearson residuals (standardised raw residuals) and the deviance residuals of our model; in both graphs, the residuals are plotted against the fitted/predicted values, which are not the fitted probabilities but the values of the linear predictor $\boldsymbol{\beta}\cdot\bf{X}$ :

A LOWESS line and a y=0 line were plotted in both graphs for easier analysis. We ideally want both LOWESS lines to be as flat and close to 0 as possible, which would mean that residuals average around 0 whatever the predicted value and there is no tendency in these residuals. We can say that the red line in the Pearson residuals graph is pretty much flat, except for the left tail leaning a little bit into the negatives, but it looks very good overall. Now, looking at the deviance residuals graph, we clearly spot a small positive slope across the entire range of predicted values. This means that there is some tendency in the deviance residuals, or that they increase on average with the predicted value. The underlying reason could be that there is a relationship between these residuals and one or more of our explanatory variables, and to find out we need to dive deeper into the behaviour of the residuals compared to each variable. But before that, let us have another look at the results above by plotting them in a different fashion:

We binned or gouped the residuals by averaging them over a certain number of predicted values. The data in these two graphs is exactly the same as in the two previous graphs, and the axes as well except that the y axis is now the average residual. We set nclass=100 which means that there are approximately 100 bins on these graphs, and given that we have 780 observations, that makes about 8 residuals averaged per bin. This view gives us another perspective and lets us decide whether we really see any tendency that would raise concern. But given the relatively small size of our data set, both of these binned residuals graphs are looking very good. For the sake of curiosity, we are still going to further investigate the deviance residuals.

We see that there is no significant trend in the graph above. We shall not forget that at both ends of the range, i.e. youngest and oldest ages, we have very few data points compared to the rest of age values closer to the centre of the range, which means that the red LOWESS line dipping into the negatives at the right end of the range is less significant. Therefore, we can consider the line as almost flat and around 0. Now let us move onto the other variables, which are all categorised variables, so we are going to need another type of graph (box plot) to visualise the behaviour of the residuals against them:

All data points were left in these box plots to better visualise where the residuals are concentrated. The solid line inside the box is the median, and it sometimes is above zero for some boxes, and below zero for the others, but what matters more is its position compared to the data points. An example of almost perfect distribution of residuals around zero is the 200m Freestyle box in the Deviance residuals vs event graph. The median is below zero, but there is only one negative residual right above it: what this means is that there are basically as many positive than negative residuals, so our model does not systematically under- or over-predict for this category. In contrast, things are looking a little bit different over in the 100m Breast_Back_Fly box: the median is negative and quite buried inside the cluster of negative residuals where it is much darker (= the biggest density of points), and this means that there are many more negative than positive residuals for this category. In other words, our model systematically over-predicts (since more residuals are negative) fasterRelay outcomes when event is 100m Breast_Back_Fly. A similar result can be observed for 100m Freestyle but the other way around: our model seems to systematically slightly under-predict fasterRelay outcomes for that category of event. We could argue that this is also happening for worlds in Deviance residuals vs meetType and the strongest swimming tradition level 3 in Deviance residuals vs swimTradition, but to a lesser extent as well.

After a quick check for interaction terms between event and other variables, the deviance residual trend did not change a bit, which unsurprisingly means that the underlying cause was not a subtle iteraction between two of the explanatory variables. However, we recall that we grouped 3 levels of event into one that now has over twice the number of observations as the two remaining levels, 100m Freestyle and 200m Freestyle. Next, we add weights to our logistic model, following that the event value of an observation is 100m Freestyle or 200m Freestyle (weight = 2) or 100m Breast_Back_Fly (weight = 1) and visualise the results:

We barely notice any change in the first box plot graph, however we notice that the red LOWESS line in Deviance residuals vs Fitted, with weights has become a little bit flatter than before (see graph to the right), especially the left tail. This is because we added more weight to observations belonging to the 100m Freestyle or 200m Freestyle categories. Indeed, we saw that the model systematically over-predicts with the 100m Breast_Back_Fly category, and after adding colours it is now clear that these values coloured in orange are distributed more to the left side of the range, which is partly responsible for the left tail of the red line dipping into the negatives. In contrast, the predicted values for 200m Freestyle coloured in dark green and 100m Freestyle coloured in purple spread further right, which is one underlying cause of the right tail of the red line going slightly up.

After adding weights to our model and colouring the residuals based on the event category of the predicted value, the reason behind the slight tendency of deviance residuals against predicted values is now clear to us: during the EDA, we saw that the proportions of fasterRelay outcomes (0/1) were different for these 3 different groups (e.g. 100m Freestyle had a much bigger percentage of 1 outcomes), and we split the character variable into 3 factors corresponding to these 3 groups being fully aware of those discrepancies – something that adding weights to the model cannot solve. Therefore, these three levels were not sampled randomly with regards to their fasterRelay outcome proportions, and the positive tendency between deviance residuals and predicted values reflects that difference in proportions. Of course this might not be the one and only reason behind that tendency, but the takeaway here is that we shouldn’t be worried about the red line not being perfectly flat because the cause behind it is now known to us and is not problematic (cannot & does not need to be fixed).

Section 2: Model validation

K-fold cross-validation

This section focuses on how well the predictive power of our model generalises to new data, or data that the model was not trained on. We also want to repeat this process several times so that we get an idea of how consistent our model is, how it performs on average, and so that our perception of its predictive power is not wronged by potential over-fitting. We are free to choose the proportions of our training and test sets, however when deciding on that we also need to take into account the number of folds we want to use and the size of dsSwim. By randomly sampling the observations of dsSwim without replacement, we shuffle the rows, and then cut the data into K folds. Given that our data set does not contain a large number of observations, we will choose a smaller K of 5. This in turn leads us to choose a 80% training/ 20% test ratio, which would mean setting one fold as the test set and the remaining data as the training set. This seems to be the best choice all things considered, since we do not have enough data to break it down into more folds, or to use 2 folds for the test set (which would leave only 3 folds for training, and that makes up too few observations for proper training).

Each fold contains 156 observations. We are now going to loop over each fold and use it as a test set, while the four other fold serve as a training set for the four prediction methods:

Let us have a quick look at what the results look like:

rsltAll is a list of 5 elements, the first being the observed values of fasterRelay for each fold, and the four others being the predicted probabilities for each method, for each fold. The next step is to extract the predicted and observed values and to summarise the various aspects of the prediction performance of the four classifier methods:

The receiver operating characteristic curve is obtained by plotting the true positive rate against the false positive rate. We can already get a quick idea of which classifier performed the best, but a much more simple way to look at it is to retrieve the area under the curve (AUC):

Statistic	Mean	St. Dev.	Min	Max
auc.Log	0.676	0.012	0.663	0.687
auc.Tree	0.610	0.009	0.601	0.626
auc.Frst	0.633	0.033	0.597	0.676
auc.Gbm	0.666	0.021	0.638	0.688

auc.All statistical summary

It now becomes obvious that the best performing classifier out of our 5-fold cross-validation is logistic regression ($AUC=0.676$), closely followed by gradient boosting ($AUC=0.666$), then by random forests ($AUC=0.633$) and finally simple classification trees ($AUC=0.610$). The classifier with the most stable and consistent predictions is classification trees ($sd=0.009$) closely followed by logistic regression ($sd=0.012$). In sum, logistic regression both performed the best on average and almost had the smallest AUC standard deviation, which means it has the best and most stable predictive performance among the four classifiers. On the ROC graph we plotted a thin dotted $y=x$ line which has an $AUC$ of 0.5: this line serves as a visual and numerical ($AUC$) reference for the predictive power of our model. If our model randomly guessed the outcome of fasterRelay, its $AUC$ would be around 0.5. The $AUC$ value obtained for logistic regression is actually not uncommon at all and can be considered acceptable, though definitely not outstanding. This once again shows that significant improvements could be made to our model to enhance our ability to predict fasterRelay, such as gathering more data (before 2015) and especially more relevant data (additional explanatory variables). Ultimately, considering that we did not carefully choose the data we collected along with solo and relay TIW, this predictive performance can definitely be considered passable enough.

$$ \sim $$

Conclusion

This project set out to find out whether elite swimmers tend to swim faster in relay events than in individual events at major swimming meets. In order to establish meaningful comparisons between these two types of races, the reaction time was subtracted from the final swim time to obtain a new variable which we baptised “time in water” (TIW): that way, individual and relay leg performances could be fairly compared with each other. We made the assumption that the resulting binary variable, which indicated whether a swimmer’s relay leg TIW was faster than their individual TIW for the same event, measured the motivation of elite swimmers to swim faster for their country’s team at major meets.

This 3rd part of the project focused on statistical modelling, where the focus was on better understanding the suspected effects discovered in our EDA of swimming event, swimmer age, national swimming tradition, type of major meet, swimmer gender, and year of the meet on faster relay TIW. These explanatory variables were not chosen on purpose, since the data for these variables was simply scraped along with individual and relay times, assuming that some kind of relationship could be established between them and faster relay TIW. For that reason, we did not really expect a model consisting of these variables to have an outstanding predictive power. All statistical tests and decisions to include or exclude variables in this part of the study were made at a 0.05 significance level, for consistency.
A first model with all the explanatory variables named above was trained on our data with logistic regression, and the goodness of fit of that model was found to be acceptable though fairly weak. It came out that the age of a swimmer at the meet and the year of the meet both had non-significant $\beta$ coefficients, which was corroborated by two $\chi^{2}$ individual effect tests with non-significant results, thus we excluded these two variables from the model. The 100m Breaststroke and 100m Butterfly factors were found to not have significantly different predictive capabilties from the 100m Backstroke factor, thus we grouped the three as a single factor. Analysing the estimated marginal probabilities not only helped in making that decision, but also revealed that the group of those three events, the 200m Freestyle, and the 100m Freestyle events all had significantly different predictive capabilities, in increasing order. We can therefore consider these three categories of swimmers as coming from different populations in terms of motivation to swim faster for their country’s team (though this conclusion has some limitations). No multicollinearity issues were found.
A second, simpler model with all the modifications described above was once again trained with logistic regression, and had a similar goodness of fit to the first model. A residual analysis was performed on this second model, and no significant concern was raised regarding the behaviour of deviance residuals: only a slight positive trend against predicted values was spotted, and the main cause was identified as being the differences in dependent variable outcomes between the three groups of events (100m Back/Breast/Fly, 200m Freestyle, 100m Freestyle) described in the previous paragraph.
Lastly, 5-fold cross-validation was performed on the second model with an 80% training / 20% holdout ratio, and logistic regression came out as the best classifier among four methods (classification trees, random forests, gradient boosting), with a passable enough ROC-AUC.

Discussion

As we said previously, there are limitations in the interpretations we made of our results related to how we performed data sampling. The small size of our data set is without a doubt one of the reasons behind the relatively disappointing goodness of fit of our model, and obtaining data from before 2015 could certainly have yielded improved results. However, given that the data for our explanatory variables was not picked following specific criteria or assumptions, the performance of our model can be considered acceptable.

Considering 100m Back/Breast/Fly, 200m Freestyle and 100m Freestyle elite swimmers as three different populations of swimmers could actually be incorrect if formulated this way, since some swimmers do swim two of these events and the corresponding relay legs (which we did not take into accountà, and would then belong to two different populations. We could instead infer that the nature of the event itself influences the difference in solo and relay TIW, and that there would be three populations of “event motivations” instead of swimmers, which is a more abstract way of looking at it. The underlying cause behind these three populations could be the nature of the relay event itself. Yannick Agnel, a former elite French freestyler, revealed that he prepared together with his 3 teammates for an entire year specifically for the 2012 London Olympics 100m Freestyle relay where they finally got the best of their American rivals. This could hint at the importance of preparing a major relay event together with the entire team to boost their motivation and increase the odds that the swimmers surpass themselves during the relay rather than in the solo event, having invested so much work and hope into preparing for that relay. In contrast, it wouldn’t be surprising if medley relay swimmers often do not train together due to the heterogenous nature of this event: two swimmers specialising in 100m Breaststroke and 100m Butterfly could be less likely to often train together, or share a similar mindset, than two 100m or 200m Freestyle swimmers. While formulating all these assumptions, the need for more data, and especially more specific data arises.

In further research, one could consult with swimming professionals to discuss potentially relevant variables to further explain the matters discussed above, such as: the mobility of the swimmer, i.e. does the swimmer often train abroad and therefore much less often with their home team? nowadays, more and more young swimmers train outside of their home country; coaching practices especially when it comes to relays; personal training habits, i.e. does the swimmer mostly train on their own with their coach, or together with other swimmers (of any nationality) – we could think that training alone reduces the probability of being more motivated by relays than solo races; etc.

Bonus: complexity control for classification trees

Based on these four plots, a minimum bucket size of 60 seems like the best choice to minimise over-fitting while maximising these four classification performance metrics.