This manuscript presents the development of novel high-order complete shape functions over star-convex polygons based on the scaled boundary finite element method. The boundary of a polygon is discretised using one-dimensional high order shape functions. Within the domain, the shape functions are analytically formulated from the equilibrium conditions of a polygon. These standard scaled boundary shape functions are augmented by introducing additional bubble functions, which renders them high-order complete up to the order of the line elements on the polygon boundary. The bubble functions are also semi-analytical and preserve the displacement compatibility between adjacent polygons. They are derived from the scaled boundary formulation by incorporating body force modes. Higher-order interpolations can be conveniently formulated by simultaneously increasing the order of the shape functions on the polygon boundary and the order of the body force mode. The resulting stiffness-matrices and mass-matrices are integrated numerically along the boundary using standard integration rules and analytically along the radial coordinate within the domain. The bubble functions improve the convergence rate of the scaled boundary finite element method in modal analyses and for problems with non-zero body forces. Numerical examples demonstrate the accuracy and convergence of the developed approach. Copyright (c) 2016 John Wiley & Sons, Ltd.