In the paper mentioned in the question, the $+U$-correction is applied to the semi-core $4d$ states of indium to increase the band gap of In2O3, which is underestimated by DFT. The band gap is important here, because the focus of the paper are charged defect formation energies, and for those the valence and conduction band edge energies are important. Otherwise, it is not too common to apply $+U$-corrections to completely filled shells.
Yes, generally, for the calculation of the formation energy of a compound, the constituents must be calculated using the same theory as the compound. In the example of the paper, bulk indium must then also be treated with LDA$+U$. There are no $U-$corrections applied to oxygen in this example, so the oxygen reference will simply be computed using LDA.
DFT is relatively good at describing delocalized metallic states, and applying $+U$-corrections to metals is problematic (the density of states is adversely affected: reduction near the Fermi level, e.g.). Therefore, as mentioned by @Tristan Maxson, at the materials project, metallic phases are computed with DFT, and transition metal oxides with partially filled $d$-shells and self-interaction-error problems with DFT$+U$. Oxide formation energies are calculated by adding a constant correction energy per $U$-corrected ion, which was fitted against a number of experimental heats of formation. This requires a constant $+U$ value for a given transition metal ion. Optimal $+U$ values will however depend on the coordination of the transition metal ion, its oxidation state, etc., and the materials project correction scheme (Jain et al., PRB 84, 045115 (2011)) has been extended to oxide-specific $U$-values (Aykol and Wolverton, PRB 90, 115105 (2014) and Voss, J. Phys. Commun. 6 035009 (2022)). With these methods one can thus compute formation energies by combining DFT and DFT$+U$ calculations, with $U$-values appropriate for a given system (even $U=0$ for metallic systems).