The production of $Upsilon(nS)$ mesons ($n=1,2,3$) in $p$Pb and Pb$p$ collisions at a centre-of-mass energy per nucleon pair $sqrt{s_{NN}}=8.16$ TeV is measured by the LHCb experiment, using a data sample corresponding to an integrated luminosity of 31.8 nb$^{-1}$. The $Upsilon(nS)$ mesons are reconstructed through their decays into two opposite-sign muons. The measurements comprise the differential production cross-sections of the $Upsilon(1S)$ and $Upsilon(2S)$ states, their forward-to-backward ratios and nuclear modification factors, performed as a function of the transverse momentum pt and rapidity in the nucleon-nucleon centre-of-mass frame $y^*$ of the $Upsilon(nS)$ states, in the kinematic range $p_{rm{T}}<25$ GeV/$c$ and $1.5<y^*<4.0$ ($-5.0<y^*<-2.5$) for $p$Pb (Pb$p$) collisions. In addition, production cross-sections for $Upsilon(3S)$ are measured integrated over phase space and the production ratios between all three $Upsilon(nS)$ states are determined. The measurements are compared to theoretical predictions and suppressions for quarkonium in $p$Pb collisions are observed.