Trace of commutators with flavor indices

First of all, your equation for gauge fields (which has the index $a$ on one side but not the other) is not correct. The closest thing that commonly comes up is begin{equation} [A^mu, A^nu] = if^{ab}_{,,,,,c} A^mu_a A^nu_b T^c end{equation} but this is a consequence of the gauge fields being valued in some Lie algebra. I.e. $A^mu = A^mu_a T^a$ where the generators satisfy begin{equation} [T^a, T^b] = if^{ab}_{,,,,,c} T^c. end{equation} The same thing is going on with the scalars in your Lagrangian. It may look like there are only as many scalar fields as values that $I$ and $J$ can take. But the partially index-free notation hides the fact that $phi^1$, $phi^2$, etc are all matrices where each entry is a scalar field.

In other words, the $phi^I$ are really $phi^I_a T^a$ meaning they are in the adjoint representation of some (global or gauge) symmetry algebra. For gauge fields, this is an absolute necessity while for scalar fields, it is a choice. And one can certainly imagine theories where the scalars are not in the adjoint. But these types of scalars (and the exact interaction you wrote) show up in one of the most famous QFTs which is $mathcal{N} = 4$ Super Yang-Mills with a $U(N)$ gauge group. In this case, $I$ and $J$ run from 1 to 6, reflecting the $SO(6)$ symmetry under which the supercharges are spinors. Holographically, this is the isometry group of the internal manifold in $AdS_5 times S^5$. But each of $phi^1, dots, phi^6$ is an $N times N$ matrix. The reason for expecting this is that they have to be in the same super multiplet (superfield) as the $A^mu$ which also are.