Probabilistic graphical models are powerful tools to quantify, visualize, and interpret network dependencies in complex biological systems such as high-throughput -omics. However, many graphical models assume sample homogeneity, limiting their effectiveness. We propose a flexible Bayesian approach called graphical regression (GraphR), which (1) incorporates sample heterogeneity at different scales through a regression-based formulation, (2) enables sparse sample-specific network estimation, (3) identifies and quantifies potential effects of heterogeneity on network structures, and (4) achieves computational efficiency via variational Bayes algorithms. We illustrate the comparative efficiency of GraphR against existing state-of-the-art methods in terms of network structure recovery and computational cost across multiple settings. We use GraphR to analyze three multi-omic and spatial transcriptomic datasets to investigate inter- and intra-sample molecular networks and delineate biological discoveries that otherwise cannot be revealed by existing approaches. We have developed a GraphR R package along with an accompanying Shiny App that provides comprehensive analysis and dynamic visualization functions.