I plot a standardized residuals against the fitted values and it does exhibit a megaphone shape. It looks like there is more variation in the lower level of fitted values.
I then conduct a Breusch-Pagan Test for heteroskedasticity.
studentized Breusch-Pagan test
data: reg4
BP = 87.88, df = 2, p-value < 2.2e-16
The p-value is quite small, which indicates that I should reject the null hypothesis and conclude heteroskedasticity. At this point, can I safely conclude I do face a heteroskedasticity and do some remedial measurements to deal with it?