Consider the case of finite dimensional Hilber space. Then may suppose that operators $T_1,\dots,T_d$ are diagonal in the same orthonormal basis $(e_1,e_2,\dots)$, denote the diagonal of $T_i$ by $(p_{i1},p_{i2},\dots)$. The square of RHS is nothing but $\sup_j |p_{1j}|^2+|p_{2j}|^2+\dots+|p_{dj}|^2$. For fixed $j$, taking $x=e_j$ we get the same thing $|p_{1j}|^2+|p_{2j}|^2+\dots+|p_{dj}|^2$ for LHS.
Now do the same in general situation. Note that
$$
\sup_{\|f\|_{L^2(\mu)}=1}\bigg[\displaystyle\sum_{i=1}^n\left(\int_Y|\varphi_i|^2|f|^2d\mu \right)\bigg]={\rm essup}\,\sum_i |\varphi_i|^2=:A.
$$
Note that (by separability of a complex plane and by the definition of essential supremum) we may find discs $D_i$ such that $\mu(\Omega)>0$, where $\Omega=\{x:\varphi_i(x)\in D_i,i=1,\dots,n\}>0$ and $\sum_i |\theta_i|^2>A-\varepsilon$ whenever $\theta_i\in D_i$. Choose any function $f$ with $\|f\|=1$ concentrated in $\Omega$. Then $\int \varphi_i |f|^2d\mu\in D_i$ (indeed, $D_i$ are convex, thus a weighted mean of elements of $D_i$ belongs to $D_i$). It implies that LHS is not less than $A-\varepsilon$.