How to approach All vs All classification problem Let's say you are building a Star trek style medical tricorder which can diagnose any medical condition. It needs to be able to detect comorbidities where a patient has multiple conditions (e.g. perhaps the patient has COVID, diabetes and lung cancer all at the same time).
How would you build a classification system to detect the most likely set of conditions?
I see two approaches:

Build one model per disease, run predictions using every model and report that the patient is suffering from any disease for which the prediction probability exceeded a threshold (e.g. 0.95)

likely_diseases = []
THRESHOLD = 0.95
for disease in disease_columns:
  model = xgb.XGBClassifier()
  model.fit(X_train, y_train[disease])
  pred_proba = model.predict_proba(patient_data)[:, 1]
  if (pred_proba &gt; THRESHOLD):
    likely_diseases.append(disease)


Build one model per combination of diseases and choose the combo with the highest probability

df['has_covid_lung_cancer_and_stroke'] = df.apply(lambda patient: patient['has_lung_cancer'] and patient['has_covid'] and patient['has_stroke'])

# create all other possible permutations of diseases

highest_probability = 0.0
most_likely_disease_combination = None
for disease_combination in disease_combinations: 
  model = xgb.XGBClassifier()
  model.fit(X_train, y_train[disease])
  pred_proba = model.predict_proba(patient_data)[:, 1]
  if (pred_proba &gt; highest_probability):
    most_likely_disease_combination = disease_combination
    highest_probability = pred_proba

It strikes me that approach two would probably be more accurate but might be so computationally expensive that it is intractable. Perhaps some pruning would take place where combinations that have exceedingly low occurrences in the training data set are discarded.
I think you can carry out a usual multiclass classification instead of manually carrying out the one-VS-all strategy (i.e. in the for loop), provided that you can generate the multilabels target(i.e. all the possible combinations of diseases you said). So, I would do something like:

generate the correct labels --&gt; if you have n possible deseases combinations, you should have in an end n target labels (so using the apply-lambda expression generates n binary columns for you, which you can also convert to a single numeric column with n possible integer values).
Below you can find a way to get your final column of target values (an integer value for each possible combination); using the wines dataset as an example where I also consider all possible combinations (here a flavour type is similar to a possible disease in your case). I guess you have a dataframe like this:


so you can create a dict mapping the possible combinations to an integer target value (for the multiclass classification):
target_labels_dict = {tuple([1, 0, 0]): 0, tuple([0, 1, 0]): 1, tuple([0, 0, 1]): 2, tuple([1, 1, 0]): 3, tuple([1, 0, 1]): 4, tuple([0, 1, 1]): 5, tuple([1, 1, 1]): 6}

and you have this:

where for instance the tuple (1, 1, 0) means 'acid' and 'sweet'...
and finally, create your target column:
target_flavours_array = np.zeros(len(y_train)).astype('int') 

for ind in y_train_df.index:
   combination_values = tuple(y_train_df.iloc[ind])
   target_flavours_array[ind] = target_labels_dict[combination_values]

getting something like:


use the softmax (or its analogous sofprob function check it out here) with the xgboost instance --&gt; below a quick example of this (training on the original 3 classes dataset for the shake of the example):

from sklearn.datasets import load_wine
from sklearn.model_selection import train_test_split
import numpy as np 
import xgboost as xgb

dataset = load_wine()
X_train, X_test, y_train, y_test = train_test_split(dataset.data, 
dataset.target, test_size=0.1, random_state=42)

dtrain = xgb.DMatrix(X_train, label=y_train)
dtest = xgb.DMatrix(X_test, label=y_test)    

params = {
    'max_depth': 4,  # maximum depth of each tree
    'eta': 0.2,  # learning rate
    'objective': 'multi:softprob',  # for multiclass problem, providing probabilities
    'num_class': len(np.unique(dataset.target))} # number of classes

iters_number = 30
xgb_model = xgb.train(params, dtrain, iters_number)
xgb_model.predict(dtest)

and, the result is the prediction probability for each possible class (diseases-combination in your use case):

here you can select the highest probability value for each row.