In finite element analysis it is well known that hexahedral elements are the preferred type of three dimensional element because of their accuracy and convergence properties. However, in general it is not possible to mesh complex geometry problems using purely hexahedral meshes. Indeed for highly complex geometries a mixture of hexahedra and tetrahedra is often required. However, in order to geometrically link hexahedra and tetrahedra, in a conforming finite element mesh, pyramid elements will be required. Until recently the basis functions of pyramid elements were not fully understood from a mathematical or computational perspective. Indeed only first-order pyramid basis functions were rigorously derived and used within the field of finite elements. This paper makes use of a method developed by Bergot that enables the generation of second and higher-order basis functions, applying them to finite element discretisations of the neutron transport equation in order to solve nuclear reactor physics, radiation shielding and nuclear criticality problems. The results demonstrate that the pyramid elements perform well in almost all cases in terms of both solution accuracy and convergence properties.