Not able properly rasterize polygon dataset using terra::rasterize()

I am attempting to create a rasterized "heat map" of tree mortality from polygon data. The polygon data consists of mapped geographies with dead trees and have an assigned value of the estimated number of trees killed. My approach to create the heat map is to first rasterize the polygon data by summing the dead tree values for all polygons within a raster pixel, and then create a percent cover raster to weight the dead tree values. The final heat map will be the product of the dead tree value raster and the percent cover raster. Finally, to smooth out the raster, I plan to downscale it to finer resolution with bilinear interpolation. However, when I follow this procedure, my final rasterized heat map does not conform to what I would expect from my polygon dataset. It seems like when I sum the dead tree values in the first rasterization, many pixels end up with a zero value that should be non-zero. Comparing it to the percent cover layer, I wonder if raster pixels that aren't completely covered end up with zero values. Is this an issue with the raster template I am using or the function I am passing to rasterize()? Do I need to change how I specify the arguments of rasterize()? Below is a reproducible example. To exaggerate the issue I made all zero values NA. In the map at the end of the example, toggle between the rasterized damage layer, the rasterized percent cover layer, and the polygon layer to see the issue.

##Loading Necessary Packages##
library(sf)
library(terra)
library(tmap)
library(tmaptools)

##For Reproducibility##
set.seed(47)

##Creating Fake Polygon Dataset##
X<-runif(250, -122.5, -121)
Y<-runif(250, 44.5, 45.5)
RAST_area<-sample(c(716.8, 204.8, 7680.0, 4096.0, 2048.0), 250, replace=TRUE)
Buff<-runif(250, 0, 10000)
DTF<-data.frame(X=X, Y=Y, RAST_area = RAST_area, Buff = Buff)
polys<-st_as_sf(DTF, coords=c(x="X", y="Y"), crs=4326) %>% 
  st_transform(6557) %>% 
  st_buffer(dist=DTF$Buff)

##Blank raster template over example data for rasterize()##
Blank<-rast(ext=st_bbox(polys), resolution=4*5280, crs="epsg:6557")

##Rasterizing the polygon data##
Poly_rast<-rasterize(x=polys, y=Blank, fun="sum", na.rm=TRUE, field="RAST_area", background=0)
Poly_rast_cov<-rasterize(x=polys, y=Blank, fun="sum", na.rm=TRUE, field="RAST_area", cover=TRUE, background=0)
Rast_4mi<-Poly_rast*Poly_rast_cov

##Resampling to finer resolution##
Rast_quarter_mi<-disagg(Rast_4mi, fact=16, method="bilinear")

##Removing Zero value pixels for mapping##
Rast_quarter_mi[Rast_quarter_mi<=0]<-NA
Poly_rast_cov[Poly_rast_cov<=0]<-NA

##Viewing the results##
tmap_mode("view")

tm_shape(polys)+
  tm_polygons(fill="RAST_area")+
tm_shape(Rast_quarter_mi)+
  tm_raster(col="RAST_area", col.scale = tm_scale_continuous_log1p(values="met.tam"))+
tm_shape(Poly_rast_cov)+
  tm_raster(col.scale=tm_scale_continuous(values="kovesi.blue_cyan"))+
tm_basemap(server=providers$OpenStreetMap)

The straightforward approach would seem to rasterize the statistic of interest at the higher resolution you desire.

You example data

library(terra)
set.seed(47)
X <- runif(250, -122.5, -121)
Y <- runif(250, 44.5, 45.5)
RAST_area <- sample(c(716.8, 204.8, 7680.0, 4096.0, 2048.0), 250, replace=TRUE)
Buff <- runif(250, 0, 10000)
DTF <- data.frame(lon=X, lat=Y, RAST_area = RAST_area, Buff = Buff)
polys <- vect(DTF, c("lon", "lat"), crs="lonlat") |> buffer(DTF$Buff) |> project("epsg:6557")

You describe your data as counts, so I think you need to first compute density

polys$area <- expanse(polys, unit="km")
polys$dead_tree_density <- polys$RAST_area / polys$area

Now rasterize

r <- rast(ext=polys, resolution=(4*5280/16), crs="epsg:6557")
r <- rasterize(x=polys, y=r, fun="sum", na.rm=TRUE, 
               field="dead_tree_density", background=NA)

plot(log(r)); lines(polys, col="gray")

If partial coverage of cells is still a concern, you could rasterize the polygons to a raster with a higher resolution, and then aggregate the results.