Abstract We study renormalization group multicritical fixed points in the ϵ-expansion of scalar field theories characterized by the symmetry of the (hyper)cubic point group HN. After reviewing the algebra of HN-invariant polynomials and arguing that there can be an entire family of multicritical (hyper)cubic solutions with ϕ2n interactions in $$ d=\frac{2n}{n-1}-\epsilon $$ d = 2 n n − 1 − ϵ dimensions, we use the general multicomponent beta functionals formalism to study the special cases d = 3 − ϵ and $$ d=\frac{8}{3}-\epsilon $$ d = 8 3 − ϵ , deriving explicitly the beta functions describing the flow of three- and four-critical (hyper)cubic models. We perform a study of their fixed points, critical exponents and quadratic deformations for various values of N, including the limit N = 0, that was reported in another paper in relation to the randomly diluted single-spin models, and an analysis of the large N limit, which turns out to be particularly interesting since it depends on the specific multicriticality. We see that, in general, the continuation in N of the random solutions is different from the continuation coming from large-N, and only the latter interpolates with the physically interesting cases of low-N such as N = 3. Finally, we also include an analysis of a theory with quintic interactions in $$ d=\frac{10}{3}-\epsilon $$ d = 10 3 − ϵ and, for completeness, the NNLO computations in d = 4 − ϵ.