diff --git a/R/anlz_fibmatrix.R b/R/anlz_fibmatrix.R
index c1dc77a0..7f4174b4 100644
--- a/R/anlz_fibmatrix.R
+++ b/R/anlz_fibmatrix.R
@@ -56,7 +56,7 @@ anlz_fibmatrix <- function(fibdata, yrrng = NULL, stas = NULL, bay_segment = NUL
                            precipdata = NULL, temporal_window = NULL, wet_threshold = NULL,
                            warn = TRUE){
 
-  geomean <- function(x){prod(x)^(1/length(x))}
+  geomean <- function(x){prod(x, na.rm = T)^(1/length(na.omit(x)))}
 
   subset_wetdry <- match.arg(subset_wetdry)
   indic <- match.arg(indic, c('fcolif', 'entero'))
@@ -207,11 +207,27 @@ anlz_fibmatrix <- function(fibdata, yrrng = NULL, stas = NULL, bay_segment = NUL
     dplyr::filter(station %in% stas) %>%
     dplyr::filter(yr >= (yrrng[1] - (lagyr - 1)) & yr <= yrrng[2]) %>%
     dplyr::filter(!is.na(indic) | indic < 0) %>%
+    tidyr::complete(
+      yr = tidyr::full_seq(c(min(yr) - (lagyr - 1), yr), 1),
+      tidyr::nesting(station)
+    ) %>%
     dplyr::summarise(
-      nyrs = length(unique(yr)),
-      .by = 'station'
+      hasdat = sum(any(!is.na(indic))),
+      .by = c('station', 'yr')
+    ) %>%
+    arrange(station, yr) %>%
+    dplyr::mutate(
+      chkyr = stats::filter(hasdat, rep(1, lagyr), sides = 1, method = 'convolution'),
+      .by = station
     ) %>%
-    dplyr::filter(nyrs >= lagyr) %>%
+    dplyr::filter(any(chkyr >= lagyr), .by = station)
+
+  # check if all stations invalid for lagyr
+  if(nrow(stasval) == 0){
+    stop('Insufficient data for lagyr')
+  }
+
+  stasval <- stasval %>%
     dplyr::pull(station) %>%
     unique()
 
@@ -224,11 +240,6 @@ anlz_fibmatrix <- function(fibdata, yrrng = NULL, stas = NULL, bay_segment = NUL
     stas <- stas[!chk]
   }
 
-  # check if all stations invalid for lagyr
-  if(sum(chk) == length(chk)){
-    stop('Insufficient data for lagyr')
-  }
-
   grp <- 'station'
   levs <- stas
   if(!is.null(bay_segment)){
@@ -243,6 +254,10 @@ anlz_fibmatrix <- function(fibdata, yrrng = NULL, stas = NULL, bay_segment = NUL
     dplyr::filter(yr >= (yrrng[1] - (lagyr - 1)) & yr <= yrrng[2]) %>%
     dplyr::filter(!is.na(indic) | indic < 0) %>%
     dplyr::rename(grp = !!grp) %>%
+    tidyr::complete(
+      yr = tidyr::full_seq(yr, 1),
+      tidyr::nesting(grp)
+    ) %>%
     summarise(
       gmean = geomean(indic),
       sumgt = sum(indic > thrsh),
diff --git a/tests/testthat/test-anlz_fibmatrix.R b/tests/testthat/test-anlz_fibmatrix.R
index 1db03dbe..6babfab4 100644
--- a/tests/testthat/test-anlz_fibmatrix.R
+++ b/tests/testthat/test-anlz_fibmatrix.R
@@ -20,13 +20,26 @@ test_that("Checking anlz_fibmatrix station error all insufficient data", {
 
 # Example data
 fibdatatst <- data.frame(
-  yr = rep(2000:2005, each = 3),
-  epchc_station = rep(letters[1:3], times = 6),
-  fcolif = runif(18, 0, 500),
-  entero = runif(18, 0, 200),
-  Latitude = runif(18, 27, 28),
-  Longitude = runif(18, -82, -81)
-)
+  epchc_station = c('a', 'b', 'c'),
+  Latitude = runif(3, 27, 28),
+  Longitude = runif(3, -82, -81),
+  area = 'FL'
+  ) %>%
+  tidyr::crossing(
+    yr = rep(2000:2005),
+    mo = rep(1:12)
+  ) %>%
+  mutate(
+    fcolif = runif(n(), 0, 500),
+    entero = runif(n(), 0, 200)
+  )
+
+test_that("Checking anlz_fibmatrix station warning for insufficient data", {
+  datchk <- fibdatatst %>%
+    filter(!(epchc_station == 'a' & yr %in% c(2000:2003)))
+  expect_warning(anlz_fibmatrix(datchk, indic = 'fcolif', stas = c('a', 'b', 'c')), regexp = 'Stations with insufficient data for lagyr: a',
+               fixed = T)
+})
 
 test_that("anlz_fibmatrix returns correct structure", {
   result <- anlz_fibmatrix(fibdatatst, indic = 'fcolif', stas = c("a", "b", "c"))