We calculate the renormalized effective two-, three- and four-body interactions for N neutral ultracold bosons in the ground state of an isotropic harmonic trap, assuming two-body interactions modeled with the combination of a zero-range and energy-dependent pseudopotential. We work to third-order in the scattering length a(t)(0) defined at zero collision energy, which is necessary to obtain both the leading-order effective four-body interaction and consistently include corrections for realistic two-body interactions. The leading-order, effective three-and four-body interaction energies are U-3 (omega) = -(0.85576...)[a(t)(0)/sigma(omega)](2) + 2.7921(1)[a(t)(0)/sigma(omega)](3) + O(a(t)(4)) and U-4(omega) = +(2.43317...)[a(t)(0)/sigma(omega)](3) + O(a(t)(4)), where omega and sigma(omega) are the harmonic oscillator frequency and length, respectively, and energies are in units of h omega. The one-standard deviation error +/- 0.0001 for the third-order coefficient in U-3(omega) is due to numerical uncertainty in estimating a slowly converging sum; the other two coefficients are either analytically or numerically exact. The effective three-and four-body interactions can play an important role in the dynamics of tightly confined and strongly correlated systems. We also performed numerical simulations for a finite-range (FR) boson-boson potential, and it was comparison to the zero-range predictions which revealed that finite-range effects must be taken into account for a realistic third-order treatment. In particular, we show that the energy-dependent pseudopotential accurately captures, through third order, the finite-range physics, and in combination with the multi-body effective interactions gives excellent agreement with the numerical simulations, validating our theoretical analysis and predictions.