How does Scikit Learn compute f1_macro for multiclass classification?

Macro-averaging doesn't work like that. A macro-average f1 score is not computed from macro-average precision and recall values.

Macro-averaging computes the value of a metric for each class and returns an unweighted average of the individual values. Thus, computing f1_score with average='macro' computes f1 scores for each class and returns the average of those scores.

If you want to compute the macro-average value yourself, specify average=None to get an array of binary f1 scores for each class, then take the mean() of that array:

binary_scores = metrics.f1_score(y_test, y_pred, average=None)
manual_f1_macro = binary_scores.mean()

Runnable demo here.

Final Update:

As a result of the very valuable comments of user2357112 (see his/her answer below too) and after reading several misunderstandings and false information on the web, finally I had to make some kind of investigation regarding macro type f1-score formula. As user2357112 below revealed also (actually first) the f1_macro has slightly different algorithm then you used at the manual calculation. Eventually I found a reliable source.

Proof that sklearn uses this formula:

A snippet from sklearn's classification.py modul's precision_recall_fscore_support() method:

    precision = _prf_divide(tp_sum, pred_sum,
                            'precision', 'predicted', average, warn_for)
    recall = _prf_divide(tp_sum, true_sum,
                         'recall', 'true', average, warn_for)
    # Don't need to warn for F: either P or R warned, or tp == 0 where pos
    # and true are nonzero, in which case, F is well-defined and zero

    f_score = ((1 + beta2) * precision * recall /
               (beta2 * precision + recall))

    f_score[tp_sum == 0] = 0.0

# Average the results

if average == 'weighted':
    weights = true_sum
    if weights.sum() == 0:
        return 0, 0, 0, None
elif average == 'samples':
    weights = sample_weight
else:
    weights = None

if average is not None:
    assert average != 'binary' or len(precision) == 1

    precision = np.average(precision, weights=weights)
    recall = np.average(recall, weights=weights)
    f_score = np.average(f_score, weights=weights)

    true_sum = None  # return no support

return precision, recall, f_score, true_sum

As we can see sklearn make the final average before precision and recall have to been averaged:

precision = np.average(precision, weights=weights)
recall = np.average(recall, weights=weights)
f_score = np.average(f_score, weights=weights)

Finally slightly altered your code:

from sklearn import datasets
from sklearn.ensemble import RandomForestClassifier
import pandas as pd
from sklearn.model_selection import train_test_split


iris = datasets.load_iris()
data=pd.DataFrame({
    'sepal length':iris.data[:,0],
    'sepal width':iris.data[:,1],
    'petal length':iris.data[:,2],
    'petal width':iris.data[:,3],
    'species':iris.target
})

X=data[['sepal length', 'sepal width', 'petal length', 'petal width']]
y=data['species']  

X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.3)
clf=RandomForestClassifier(n_estimators=100)

clf.fit(X_train,y_train)
y_pred=clf.predict(X_test)

#Compute metrics using scikit
from sklearn import metrics
print(metrics.confusion_matrix(y_test, y_pred))
print(metrics.classification_report(y_test, y_pred))
pre_macro = metrics.precision_score(y_test, y_pred, average="macro")
recall_macro = metrics.recall_score(y_test, y_pred, average="macro")
f1_macro_scikit = metrics.f1_score(y_test, y_pred, average="macro")

f1_score_raw = metrics.f1_score(y_test, y_pred, average=None)

f1_macro_manual = f1_score_raw.mean()

print ("Prec_macro_scikit:", pre_macro)
print ("Rec_macro_scikit:", recall_macro)
print ("f1_macro_scikit:", f1_macro_scikit)

print("f1_score_raw:", f1_score_raw)
print("f1_macro_manual:", f1_macro_manual)

Out:

[[16  0  0]
 [ 0 15  0]
 [ 0  6  8]]
             precision    recall  f1-score   support

          0       1.00      1.00      1.00        16
          1       0.71      1.00      0.83        15
          2       1.00      0.57      0.73        14

avg / total       0.90      0.87      0.86        45

Prec_macro_scikit: 0.9047619047619048
Rec_macro_scikit: 0.8571428571428571
f1_macro_scikit: 0.8535353535353535
f1_score_raw: [1.         0.83333333 0.72727273]
f1_macro_manual: 0.8535353535353535

Alternatively you can make a "manual calculation" similarly like you did:

import numpy as np

pre = metrics.precision_score(y_test, y_pred, average=None)
recall = metrics.recall_score(y_test, y_pred, average=None)

f1_macro_manual = 2 * pre * recall / (pre + recall )
f1_macro_manual = np.average(f1_macro_manual)

print("f1_macro_manual_2:", f1_macro_manual)

Out:

f1_macro_manual_2: 0.8535353535353535