Graphene nanoclusters (GNCs), part of the graphene family of nanomaterials, are attracting interest as promising Pt free catalysts for the Oxygen Reduction Reaction (ORR) at the cathode of the fuel cells. This investigation focuses on the theoretical prediction of their possible edge activity towards ORR determined by their shape, size and termination type. In this framework, using density functional theory (DFT) we examined the zigzag (zz) and armchair terminated (ac) triangular (T), rhombohedral (R) and hexagonal (H) GNCs of various sizes, ranging from C13 to C114. The theoretical onset overpotential (ηth,onset) indicate that several GNC edges could be active for ORR. The model indicates that the zigzag edges of all triangular GNCs (C13–C46) as well as of the largest investigated rhombohedral shapes (C48, C70) are active for ORR. A correlation between the band gaps, spin densities and p-band centers of the edge carbon atoms and the adsorption energies of the ORR intermediate moieties (HOO*/O*/HO*) are discussed and which are directly related to their activity. One trend is related to the decrease of the adsorption energies with the decrease of the band gap and accordingly an increase of their ORR activities. The case studies, show that the explicit water model developed by using machine learning potentials and combined with DFT, stabilizes the ORR intermediates stronger than the implicit one, thus increasing slightly the onset overpotential. Overall, this computational study enhances our fundamental understanding of the edges of GNC for their potential application in electrocatalysis without doping.