diff --git a/Tests/test_imagestat.py b/Tests/test_imagestat.py
index 0baab7ce21e..47e534696aa 100644
--- a/Tests/test_imagestat.py
+++ b/Tests/test_imagestat.py
@@ -67,3 +67,9 @@ def test_zero_count() -> None:
     assert st.mean == [0]
     assert st.rms == [0]
     assert st.var == [0]
+
+
+def test_variance_rounding() -> None:
+    im = Image.new("L", (919, 405), 255)
+    st = ImageStat.Stat(im)
+    assert st.var == [0.0]
diff --git a/docs/releasenotes/13.0.0.rst b/docs/releasenotes/13.0.0.rst
index 46dbed3b58b..1f2ee494f01 100644
--- a/docs/releasenotes/13.0.0.rst
+++ b/docs/releasenotes/13.0.0.rst
@@ -139,3 +139,9 @@ Fixed ImageChops.offset() for 16-bit images
 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
 
 :py:meth:`~PIL.ImageChops.offset` now works correctly for 16-bit-per-channel images.
+
+Fixed ImageStat variance rounding
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+Fixed a floating-point rounding error that may have produced a negative variance in
+:py:attr:`~PIL.ImageStat.Stat.var`.
diff --git a/src/PIL/ImageStat.py b/src/PIL/ImageStat.py
index fde0256b2f8..769c035b79b 100644
--- a/src/PIL/ImageStat.py
+++ b/src/PIL/ImageStat.py
@@ -96,26 +96,26 @@ def count(self) -> list[int]:
         return [sum(self.h[i : i + 256]) for i in range(0, len(self.h), 256)]
 
     @cached_property
-    def sum(self) -> list[float]:
+    def sum(self) -> list[int]:
         """Sum of all pixels for each band in the image."""
 
         v = []
         for i in range(0, len(self.h), 256):
-            layer_sum = 0.0
+            layer_sum = 0
             for j in range(256):
                 layer_sum += j * self.h[i + j]
             v.append(layer_sum)
         return v
 
     @cached_property
-    def sum2(self) -> list[float]:
+    def sum2(self) -> list[int]:
         """Squared sum of all pixels for each band in the image."""
 
         v = []
         for i in range(0, len(self.h), 256):
-            sum2 = 0.0
+            sum2 = 0
             for j in range(256):
-                sum2 += (j**2) * float(self.h[i + j])
+                sum2 += (j**2) * self.h[i + j]
             v.append(sum2)
         return v
 
@@ -153,7 +153,7 @@ def var(self) -> list[float]:
         """Variance for each band in the image."""
         return [
             (
-                (self.sum2[i] - (self.sum[i] ** 2.0) / self.count[i]) / self.count[i]
+                (self.sum2[i] - (self.sum[i] ** 2) / self.count[i]) / self.count[i]
                 if self.count[i]
                 else 0
             )