diff --git a/DESCRIPTION b/DESCRIPTION
index 8af97ee..1c550ef 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -25,7 +25,9 @@ Imports:
     jsonlite,
     nhdplusTools,
     sf,
-    methods
+    methods,
+    tidygraph,
+    igraph
 Suggests:
     testthat
 License: GPL (>= 3)
@@ -39,6 +41,7 @@ Collate:
     'OGRSQLDriver.R'
     'OGRSQLResult.R'
     'find_origin.R'
+    'get_shortest_path.R'
     'get_subset.R'
     'hfsubsetR-package.R'
     'query.R'
diff --git a/NAMESPACE b/NAMESPACE
index 2cd8db1..cffea8b 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -9,6 +9,7 @@ export(as_ogr)
 export(collect)
 export(find_origin)
 export(get_hydrofabric)
+export(get_shortest_path)
 export(get_subset)
 export(get_vpu_fabric)
 export(query)
@@ -57,6 +58,7 @@ importFrom(DBI,dbUnloadDriver)
 importFrom(arrow,open_dataset)
 importFrom(dbplyr,dbplyr_edition)
 importFrom(dbplyr,sql)
+importFrom(dplyr,"%>%")
 importFrom(dplyr,`%>%`)
 importFrom(dplyr,any_of)
 importFrom(dplyr,collect)
@@ -73,6 +75,10 @@ importFrom(glue,glue)
 importFrom(httr,GET)
 importFrom(httr,progress)
 importFrom(httr,write_disk)
+importFrom(igraph,E)
+importFrom(igraph,V)
+importFrom(igraph,as.igraph)
+importFrom(igraph,shortest_paths)
 importFrom(jsonlite,toJSON)
 importFrom(methods,callNextMethod)
 importFrom(methods,new)
@@ -89,7 +95,11 @@ importFrom(sf,st_bbox)
 importFrom(sf,st_crs)
 importFrom(sf,st_layers)
 importFrom(sf,st_point)
+importFrom(sf,st_read)
 importFrom(sf,st_set_crs)
 importFrom(sf,st_sfc)
+importFrom(sf,st_write)
 importFrom(sf,write_sf)
+importFrom(tidygraph,activate)
+importFrom(tidygraph,as_tbl_graph)
 importFrom(utils,packageVersion)
diff --git a/R/get_shortest_path.R b/R/get_shortest_path.R
new file mode 100644
index 0000000..a12aed5
--- /dev/null
+++ b/R/get_shortest_path.R
@@ -0,0 +1,169 @@
+#' @title Find the shortest path between two points in a hydrofabric network
+#' @param start_id the starting NHDPlusV2 COMID id of the shortest path. datatype: int / vector of int e.g., 61297116 or c(61297116 , 6129261) 
+#' @param end_id the ending NHDPlusV2 COMID id of the shortest path. datatype: int / vector of int e.g., 61297116 or c(61297116 , 6129261) 
+#' @param gpkg a local gpkg file
+#' @param filename If filename is provided, data will be written using the filename
+#' @param lyrs layers to extract
+#' @return An sf object containing the shortest path
+#' @details This function identifies the shortest path between two nodes in a hydrofabric network.
+#' @author Caleb Novinger  <Caleb.D.Novinger@rtx.com>
+#' @author Tadd Bindas <Tadd.N.Bindas@rtx.com>
+#' @export
+#' @importFrom sf st_read st_write
+#' @importFrom tidygraph as_tbl_graph activate
+#' @importFrom dplyr filter mutate pull %>%
+#' @importFrom igraph as.igraph V E shortest_paths
+#' 
+
+get_shortest_path <- function(
+    start_id,
+    end_id,
+    gpkg = NULL,
+    filename = NULL,
+    lyrs = c("divides", "flowpaths", "network", "nexus", "flowpath-attributes")
+  ) {
+  
+    local_subset <- get_subset(
+      comid=end_id, 
+      gpkg=gpkg,
+      lyrs=lyrs
+    )
+    
+    network <- local_subset$network
+
+    start_node <- network %>% 
+        filter(hf_id == start_id) %>% 
+        pull(id) %>%
+        .[1] 
+    end_node <- network %>% 
+        filter(hf_id == end_id) %>% 
+        pull(id) %>%
+        .[1] 
+    
+    weights <- network$lengthkm
+    weights[is.na(weights)] <- 1
+    
+    topo_graph <- igraph::graph_from_data_frame(
+      d = data.frame(from = network$id, to = network$toid, weight = weights),
+      directed = TRUE
+    )
+    
+    all_nodes <- igraph::V(topo_graph)$name  
+    if (!(start_node %in% all_nodes)) {
+      stop(paste("ERROR: Start node", start_id, "not found in the graph!"))
+    }
+    if (!(end_node %in% all_nodes)) {
+      stop(paste("ERROR: End node", end_id, "not found in the graph!"))
+    }
+    
+    shortest_path <- igraph::shortest_paths(
+      topo_graph,
+      from = which(igraph::V(topo_graph)$name == start_node),
+      to = which(igraph::V(topo_graph)$name == end_node),
+      weights = igraph::E(topo_graph)$weight,
+      mode = "out"
+    )
+    
+    path_nodes <- igraph::V(topo_graph)[shortest_path$vpath[[1]]]$name
+    
+    result <- list()
+    
+    if ("flowpaths" %in% lyrs || length(lyrs) == 0) {
+      flowpaths <- local_subset$flowpaths
+      flowpaths$id <- as.character(flowpaths$id)
+      shortest_flowpaths <- flowpaths[flowpaths$id %in% path_nodes, ]
+      result$flowpaths <- shortest_flowpaths
+    }
+    
+    if ("divides" %in% lyrs || length(lyrs) == 0) {
+      if (exists("shortest_flowpaths", inherits = FALSE)) {
+        divides <- local_subset$divides
+        divides$id <- as.character(divides$id)
+        shortest_divides <- divides[divides$id %in% shortest_flowpaths$id, ]
+        result$divides <- shortest_divides
+      }
+    }
+    
+    if ("nexus" %in% lyrs || length(lyrs) == 0) {
+      if (exists("shortest_flowpaths", inherits = FALSE)) {
+        nexus <- local_subset$nexus
+        nexus$id <- as.character(nexus$id)
+        shortest_nexus <- nexus[nexus$id %in% shortest_flowpaths$toid, ]
+        result$nexus <- shortest_nexus
+      }
+    }
+    
+    if ("network" %in% lyrs || length(lyrs) == 0) {
+      network <- local_subset$network
+      network$id <- as.character(network$id)
+      network$toid <- as.character(network$toid)
+      shortest_network <- network[network$id %in% path_nodes | 
+                                    network$toid %in% path_nodes, ]
+      result$network <- shortest_network
+    }
+    
+    if ("flowpath-attributes" %in% lyrs || length(lyrs) == 0) {
+      if (exists("shortest_flowpaths", inherits = FALSE)) {
+        flowpath_attributes <- local_subset$`flowpath-attributes`
+        flowpath_attributes$id <- as.character(flowpath_attributes$id)
+        shortest_flowpath_attributes <- flowpath_attributes[flowpath_attributes$id %in% shortest_flowpaths$id, ]
+        result$`flowpath-attributes` <- shortest_flowpath_attributes
+      }
+    }
+    
+    if ("flowpath-attributes-ml" %in% lyrs || length(lyrs) == 0) {
+      if (exists("shortest_flowpaths", inherits = FALSE)) {
+        flowpath_attributes_ml <- local_subset$`flowpath-attributes-ml`
+        flowpath_attributes_ml$id <- as.character(flowpath_attributes_ml$id)
+        shortest_flowpath_attributes_ml <- flowpath_attributes_ml[flowpath_attributes_ml$id %in% shortest_flowpaths$id, ]
+        result$`flowpath-attributes-ml` <- shortest_flowpath_attributes_ml
+      }
+    }
+    
+    if ("hydrolocations" %in% lyrs || length(lyrs) == 0) {
+      if (exists("shortest_flowpaths", inherits = FALSE)) {
+        hydrolocations <- local_subset$hydrolocations
+        hydrolocations$id <- as.character(hydrolocations$id)
+        shortest_hydrolocations <- hydrolocations[hydrolocations$id %in% shortest_flowpaths$id, ]
+        result$hydrolocations <- shortest_hydrolocations
+      }
+    }
+    
+    if ("lakes" %in% lyrs || length(lyrs) == 0) {
+      if (exists("network", inherits = FALSE)) {
+        lakes <- local_subset$lakes
+        shortest_lakes <- lakes[lakes$hf_id %in% network$hf_id, ]
+        result$lakes <- shortest_lakes
+      }
+    }
+    
+    if ("pois" %in% lyrs || length(lyrs) == 0) {
+      if (exists("shortest_flowpaths", inherits = FALSE)) {
+        pois <- local_subset$pois
+        pois$id <- as.character(pois$id)
+        shortest_pois <- pois[pois$id %in% shortest_flowpaths$id, ]
+        result$pois <- shortest_pois
+      }
+    }
+    
+    if ("divide-attributes" %in% lyrs || length(lyrs) == 0) {
+      if (exists("shortest_divides", inherits = FALSE)) {
+        divide_attributes <- local_subset$`divide-attributes`
+        divide_attributes$divide_id <- as.character(divide_attributes$divide_id)
+        shortest_divide_attributes <- divide_attributes[divide_attributes$divide_id %in% shortest_divides$divide_id, ]
+        result$`divide-attributes` <- shortest_divide_attributes
+      }
+    }
+    
+    if (!is.null(filename)) {
+      gpkg_dir <- dirname(gpkg)
+      output_filename <- paste0(filename, ".gpkg")
+      output_gpkg_path <- file.path(gpkg_dir, output_filename)
+      
+      for (layer_name in names(result)) {
+        sf::st_write(result[[layer_name]], output_gpkg_path, layer = layer_name, append = FALSE)
+      }
+    }
+    
+    return(result)
+}
\ No newline at end of file
diff --git a/man/get_shortest_path.Rd b/man/get_shortest_path.Rd
new file mode 100644
index 0000000..5a504a3
--- /dev/null
+++ b/man/get_shortest_path.Rd
@@ -0,0 +1,39 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/get_shortest_path.R
+\name{get_shortest_path}
+\alias{get_shortest_path}
+\title{Find the shortest path between two points in a hydrofabric network}
+\usage{
+get_shortest_path(
+  start_id,
+  end_id,
+  gpkg = NULL,
+  filename = NULL,
+  lyrs = c("divides", "flowpaths", "network", "nexus", "flowpath-attributes")
+)
+}
+\arguments{
+\item{start_id}{the starting NHDPlusV2 COMID id of the shortest path. datatype: int / vector of int e.g., 61297116 or c(61297116 , 6129261)}
+
+\item{end_id}{the ending NHDPlusV2 COMID id of the shortest path. datatype: int / vector of int e.g., 61297116 or c(61297116 , 6129261)}
+
+\item{gpkg}{a local gpkg file}
+
+\item{filename}{If filename is provided, data will be written using the filename}
+
+\item{lyrs}{layers to extract}
+}
+\value{
+An sf object containing the shortest path
+}
+\description{
+Find the shortest path between two points in a hydrofabric network
+}
+\details{
+This function identifies the shortest path between two nodes in a hydrofabric network.
+}
+\author{
+Caleb Novinger  <Caleb.D.Novinger@rtx.com>
+
+Tadd Bindas <Tadd.N.Bindas@rtx.com>
+}