In this work we study vacuum decay and bubble nucleation in models of $f(R)$ higher curvature gravity. Building upon the analysis of Coleman-De Luccia (CDL), we present the formalism to calculate the Euclidean action and the bounce solution for a general $f(R)$ gravity in the thin wall approximation. We calculate the size of the nucleated bubble and the decay exponent for the Starobinsky model and its higher power extensions. We have shown that in the Starobinsky model with a typical potential the nucleated bubble has a larger size in comparison to the CDL bubble and with a lower tunneling rate. However, for higher power extension of the Starobinsky model the size of the bubble and the tunneling exponent can be larger or smaller than the CDL bubble depending on the model parameters. As a counterintuitive example, we have shown that a bubble with a larger size than the CDL bubble but with a higher nucleation rate can be formed in $f(R)$ gravity.