expression_filepath <- "fcortex.expression.txt"
coordinates_filepath <- "fcortex.coordinates.txt"
genes_filepath <- "genes"

expr <- readExpression(expression_filepath)
coords <- readCoordinates(coordinates_filepath)
fields <- coords$field
coords <- coords$coords
genes <- readGenes(genes_filepath)

spatial_network <- createSpatialNetwork(coords = coords, k = 4)

hmrfParams <- initializeHMRF(expr = expr, 
                             spatial_network = spatial_network, 
                             spatial_genes = genes, 
                             k = 3)

hmrfResults <- runHMRF(expr = expr, 
                       nei = hmrfParams$nei, 
                       numnei = hmrfParams$numnei, 
                       blocks = hmrfParams$blocks, 
                       damp = hmrfParams$damp, 
                       mu = hmrfParams$mu, 
                       sigma = hmrfParams$sigma, 
                       edgelist = hmrfParams$edgelist, 
                       k = 3, 
                       genes = genes, 
                       beta_initial_value = 0, 
                       beta_increment_value = 10, 
                       beta_iterations = 5)

hmrfPlots <- plotHMRF(spatial_network = spatial_network, HMRFoutput = hmrfResults)
hmrfPlots