diff --git a/Tests/test_image_filter.py b/Tests/test_image_filter.py
index 5a22738a9..e566cd055 100644
--- a/Tests/test_image_filter.py
+++ b/Tests/test_image_filter.py
@@ -35,12 +35,14 @@ from .helper import assert_image_equal, hopper
         ImageFilter.UnsharpMask(10),
     ),
 )
-@pytest.mark.parametrize("mode", ("L", "I", "RGB", "CMYK"))
+@pytest.mark.parametrize(
+    "mode", ("L", "I", "I;16", "I;16L", "I;16B", "I;16N", "RGB", "CMYK")
+)
 def test_sanity(
     filter_to_apply: ImageFilter.Filter | type[ImageFilter.Filter], mode: str
 ) -> None:
     im = hopper(mode)
-    if mode != "I" or (
+    if mode[0] != "I" or (
         callable(filter_to_apply)
         and issubclass(filter_to_apply, ImageFilter.BuiltinFilter)
     ):
@@ -49,7 +51,9 @@ def test_sanity(
         assert out.size == im.size
 
 
-@pytest.mark.parametrize("mode", ("L", "I", "RGB", "CMYK"))
+@pytest.mark.parametrize(
+    "mode", ("L", "I", "I;16", "I;16L", "I;16B", "I;16N", "RGB", "CMYK")
+)
 def test_sanity_error(mode: str) -> None:
     im = hopper(mode)
     with pytest.raises(TypeError):
@@ -150,7 +154,9 @@ def test_kernel_not_enough_coefficients() -> None:
         ImageFilter.Kernel((3, 3), (0, 0))
 
 
-@pytest.mark.parametrize("mode", ("L", "LA", "I", "RGB", "CMYK"))
+@pytest.mark.parametrize(
+    "mode", ("L", "LA", "I", "I;16", "I;16L", "I;16B", "I;16N", "RGB", "CMYK")
+)
 def test_consistency_3x3(mode: str) -> None:
     with Image.open("Tests/images/hopper.bmp") as source:
         with Image.open("Tests/images/hopper_emboss.bmp") as reference:
@@ -166,7 +172,9 @@ def test_consistency_3x3(mode: str) -> None:
             assert_image_equal(source.filter(kernel), reference)
 
 
-@pytest.mark.parametrize("mode", ("L", "LA", "I", "RGB", "CMYK"))
+@pytest.mark.parametrize(
+    "mode", ("L", "LA", "I", "I;16", "I;16L", "I;16B", "I;16N", "RGB", "CMYK")
+)
 def test_consistency_5x5(mode: str) -> None:
     with Image.open("Tests/images/hopper.bmp") as source:
         with Image.open("Tests/images/hopper_emboss_more.bmp") as reference:
diff --git a/src/libImaging/Filter.c b/src/libImaging/Filter.c
index fbd6b425f..7b7b2e429 100644
--- a/src/libImaging/Filter.c
+++ b/src/libImaging/Filter.c
@@ -26,6 +26,8 @@
 
 #include "Imaging.h"
 
+#define ROUND_UP(f) ((int)((f) >= 0.0 ? (f) + 0.5F : (f) - 0.5F))
+
 static inline UINT8
 clip8(float in) {
     if (in <= 0.0) {
@@ -105,6 +107,22 @@ ImagingExpand(Imaging imIn, int xmargin, int ymargin) {
     return imOut;
 }
 
+float
+kernel_i16(int size, UINT8 *in0, int x, const float *kernel, int bigendian) {
+    int i;
+    float result = 0;
+    int half_size = (size - 1) / 2;
+    for (i = 0; i < size; i++) {
+        int x1 = x + i - half_size;
+        result += _i2f(
+                      in0[x1 * 2 + (bigendian ? 1 : 0)] +
+                      (in0[x1 * 2 + (bigendian ? 0 : 1)] >> 8)
+                  ) *
+                  kernel[i];
+    }
+    return result;
+}
+
 void
 ImagingFilter3x3(Imaging imOut, Imaging im, const float *kernel, float offset) {
 #define KERNEL1x3(in0, x, kernel, d)                               \
@@ -135,6 +153,16 @@ ImagingFilter3x3(Imaging imOut, Imaging im, const float *kernel, float offset) {
                 out[x] = in0[x];
             }
         } else {
+            int bigendian = 0;
+            if (im->type == IMAGING_TYPE_SPECIAL) {
+                if (strcmp(im->mode, "I;16B") == 0
+#ifdef WORDS_BIGENDIAN
+                    || strcmp(im->mode, "I;16N") == 0
+#endif
+                ) {
+                    bigendian = 1;
+                }
+            }
             for (y = 1; y < im->ysize - 1; y++) {
                 UINT8 *in_1 = (UINT8 *)im->image[y - 1];
                 UINT8 *in0 = (UINT8 *)im->image[y];
@@ -142,14 +170,31 @@ ImagingFilter3x3(Imaging imOut, Imaging im, const float *kernel, float offset) {
                 UINT8 *out = (UINT8 *)imOut->image[y];
 
                 out[0] = in0[0];
+                if (im->type == IMAGING_TYPE_SPECIAL) {
+                    out[1] = in0[1];
+                }
                 for (x = 1; x < im->xsize - 1; x++) {
                     float ss = offset;
-                    ss += KERNEL1x3(in1, x, &kernel[0], 1);
-                    ss += KERNEL1x3(in0, x, &kernel[3], 1);
-                    ss += KERNEL1x3(in_1, x, &kernel[6], 1);
-                    out[x] = clip8(ss);
+                    if (im->type == IMAGING_TYPE_SPECIAL) {
+                        ss += kernel_i16(3, in1, x, &kernel[0], bigendian);
+                        ss += kernel_i16(3, in0, x, &kernel[3], bigendian);
+                        ss += kernel_i16(3, in_1, x, &kernel[6], bigendian);
+                        int ss_int = ROUND_UP(ss);
+                        out[x * 2 + (bigendian ? 1 : 0)] = clip8(ss_int % 256);
+                        out[x * 2 + (bigendian ? 0 : 1)] = clip8(ss_int >> 8);
+                    } else {
+                        ss += KERNEL1x3(in1, x, &kernel[0], 1);
+                        ss += KERNEL1x3(in0, x, &kernel[3], 1);
+                        ss += KERNEL1x3(in_1, x, &kernel[6], 1);
+                        out[x] = clip8(ss);
+                    }
+                }
+                if (im->type == IMAGING_TYPE_SPECIAL) {
+                    out[x * 2] = in0[x * 2];
+                    out[x * 2 + 1] = in0[x * 2 + 1];
+                } else {
+                    out[x] = in0[x];
                 }
-                out[x] = in0[x];
             }
         }
     } else {
@@ -261,6 +306,16 @@ ImagingFilter5x5(Imaging imOut, Imaging im, const float *kernel, float offset) {
                 out[x + 1] = in0[x + 1];
             }
         } else {
+            int bigendian = 0;
+            if (im->type == IMAGING_TYPE_SPECIAL) {
+                if (strcmp(im->mode, "I;16B") == 0
+#ifdef WORDS_BIGENDIAN
+                    || strcmp(im->mode, "I;16N") == 0
+#endif
+                ) {
+                    bigendian = 1;
+                }
+            }
             for (y = 2; y < im->ysize - 2; y++) {
                 UINT8 *in_2 = (UINT8 *)im->image[y - 2];
                 UINT8 *in_1 = (UINT8 *)im->image[y - 1];
@@ -271,17 +326,39 @@ ImagingFilter5x5(Imaging imOut, Imaging im, const float *kernel, float offset) {
 
                 out[0] = in0[0];
                 out[1] = in0[1];
+                if (im->type == IMAGING_TYPE_SPECIAL) {
+                    out[2] = in0[2];
+                    out[3] = in0[3];
+                }
                 for (x = 2; x < im->xsize - 2; x++) {
                     float ss = offset;
-                    ss += KERNEL1x5(in2, x, &kernel[0], 1);
-                    ss += KERNEL1x5(in1, x, &kernel[5], 1);
-                    ss += KERNEL1x5(in0, x, &kernel[10], 1);
-                    ss += KERNEL1x5(in_1, x, &kernel[15], 1);
-                    ss += KERNEL1x5(in_2, x, &kernel[20], 1);
-                    out[x] = clip8(ss);
+                    if (im->type == IMAGING_TYPE_SPECIAL) {
+                        ss += kernel_i16(5, in2, x, &kernel[0], bigendian);
+                        ss += kernel_i16(5, in1, x, &kernel[5], bigendian);
+                        ss += kernel_i16(5, in0, x, &kernel[10], bigendian);
+                        ss += kernel_i16(5, in_1, x, &kernel[15], bigendian);
+                        ss += kernel_i16(5, in_2, x, &kernel[20], bigendian);
+                        int ss_int = ROUND_UP(ss);
+                        out[x * 2 + (bigendian ? 1 : 0)] = clip8(ss_int % 256);
+                        out[x * 2 + (bigendian ? 0 : 1)] = clip8(ss_int >> 8);
+                    } else {
+                        ss += KERNEL1x5(in2, x, &kernel[0], 1);
+                        ss += KERNEL1x5(in1, x, &kernel[5], 1);
+                        ss += KERNEL1x5(in0, x, &kernel[10], 1);
+                        ss += KERNEL1x5(in_1, x, &kernel[15], 1);
+                        ss += KERNEL1x5(in_2, x, &kernel[20], 1);
+                        out[x] = clip8(ss);
+                    }
+                }
+                if (im->type == IMAGING_TYPE_SPECIAL) {
+                    out[x * 2 + 0] = in0[x * 2 + 0];
+                    out[x * 2 + 1] = in0[x * 2 + 1];
+                    out[x * 2 + 2] = in0[x * 2 + 2];
+                    out[x * 2 + 3] = in0[x * 2 + 3];
+                } else {
+                    out[x + 0] = in0[x + 0];
+                    out[x + 1] = in0[x + 1];
                 }
-                out[x + 0] = in0[x + 0];
-                out[x + 1] = in0[x + 1];
             }
         }
     } else {
@@ -383,7 +460,8 @@ ImagingFilter(Imaging im, int xsize, int ysize, const FLOAT32 *kernel, FLOAT32 o
     Imaging imOut;
     ImagingSectionCookie cookie;
 
-    if (im->type != IMAGING_TYPE_UINT8 && im->type != IMAGING_TYPE_INT32) {
+    if (im->type == IMAGING_TYPE_FLOAT32 ||
+        (im->type == IMAGING_TYPE_SPECIAL && im->bands != 1)) {
         return (Imaging)ImagingError_ModeError();
     }