Weighted Sweep Example Script¶

Here, we sweep over the alpha parameter, which is the value that determines communicability in a weighted GNM

In [2]:

# basic imports from within the package plus torch and time
from gnm import defaults, utils, fitting, generative_rules, weight_criteria, evaluation
import torch
import time
import numpy as np

# Remember, the device is the hardware on which the code will run. Use a GPU for speed if 
# you have one, as this utlizes parallelization (running lots of models at the same time).
DEVICE = torch.device("cuda:0" if torch.cuda.is_available() else "cpu")

# load the example data
distance_matrix = defaults.get_distance_matrix(device=DEVICE)
weighted_consensus_network = defaults.get_weighted_network(device=DEVICE)
binary_consensus_network = defaults.get_binary_network(device=DEVICE)

# basic imports from within the package plus torch and time from gnm import defaults, utils, fitting, generative_rules, weight_criteria, evaluation import torch import time import numpy as np # Remember, the device is the hardware on which the code will run. Use a GPU for speed if # you have one, as this utlizes parallelization (running lots of models at the same time). DEVICE = torch.device("cuda:0" if torch.cuda.is_available() else "cpu") # load the example data distance_matrix = defaults.get_distance_matrix(device=DEVICE) weighted_consensus_network = defaults.get_weighted_network(device=DEVICE) binary_consensus_network = defaults.get_binary_network(device=DEVICE)

Parameter setting + Defining the Sweep¶

In [ ]:

# set fixed eta and gamma - usually you'd want to sweep through these too but here we will keep them
# constant for ease of testing
eta = torch.Tensor([-0.1])
gamma = torch.Tensor([0.1])

# set space for weighted sweep - test for 2x2x2 grid search,
# uncomment np.linspace for a more extensive search
alpha_values = [1, 2, 3] # np.linspace(-1, 1, 0.1)
omega_values = [1, 2, 3] # np.linspace(-1, 1, 0.1)

# define the number of simulations we want to run
num_simulations = 100

# set the number of connections - remember to count per connection without weight initially
num_connections = int( torch.where(weighted_consensus_network > 1, 1, 0).sum().item() / 2 )
print(f"The weighted consensus network contains {num_connections} connections.")

# set fixed eta and gamma - usually you'd want to sweep through these too but here we will keep them # constant for ease of testing eta = torch.Tensor([-0.1]) gamma = torch.Tensor([0.1]) # set space for weighted sweep - test for 2x2x2 grid search, # uncomment np.linspace for a more extensive search alpha_values = [1, 2, 3] # np.linspace(-1, 1, 0.1) omega_values = [1, 2, 3] # np.linspace(-1, 1, 0.1) # define the number of simulations we want to run num_simulations = 100 # set the number of connections - remember to count per connection without weight initially num_connections = int( torch.where(weighted_consensus_network > 1, 1, 0).sum().item() / 2 ) print(f"The weighted consensus network contains {num_connections} connections.")

The weighted consensus network contains 400 connections.

In [ ]:

# Even though we're running a weighted model, this operates on top of a binary network so we
# must define the binary sweep parameters too. 
binary_sweep_parameters = fitting.BinarySweepParameters(
    eta = eta,
    gamma = gamma,
    lambdah = torch.Tensor([0.0]),
    distance_relationship_type = ["powerlaw"],
    preferential_relationship_type = ["powerlaw"],
    heterochronicity_relationship_type = ["powerlaw"],
    generative_rule = [generative_rules.MatchingIndex()],
    num_iterations = [num_connections],
)

# we create a set of optimization criteria using our many alpha values - the question
# we're asking here is 'What is the optimium level of communicability that best matches the 
# given real life connectome?'
weighted_sweep_parameters = fitting.WeightedSweepParameters(
    alpha = alpha_values,
    optimisation_criterion = [ 
        weight_criteria.DistanceWeightedCommunicability(distance_matrix=distance_matrix, omega=omega_value)
        for omega_value in omega_values],
)

# Now we define the sweep configuation - this allows for running through every single combination 
# of parameters and saving the results.
sweep_config = fitting.SweepConfig(
    binary_sweep_parameters = binary_sweep_parameters,
    weighted_sweep_parameters = weighted_sweep_parameters,
    num_simulations = num_simulations,
    distance_matrix = [distance_matrix],
)

# Even though we're running a weighted model, this operates on top of a binary network so we # must define the binary sweep parameters too. binary_sweep_parameters = fitting.BinarySweepParameters( eta = eta, gamma = gamma, lambdah = torch.Tensor([0.0]), distance_relationship_type = ["powerlaw"], preferential_relationship_type = ["powerlaw"], heterochronicity_relationship_type = ["powerlaw"], generative_rule = [generative_rules.MatchingIndex()], num_iterations = [num_connections], ) # we create a set of optimization criteria using our many alpha values - the question # we're asking here is 'What is the optimium level of communicability that best matches the # given real life connectome?' weighted_sweep_parameters = fitting.WeightedSweepParameters( alpha = alpha_values, optimisation_criterion = [ weight_criteria.DistanceWeightedCommunicability(distance_matrix=distance_matrix, omega=omega_value) for omega_value in omega_values], ) # Now we define the sweep configuation - this allows for running through every single combination # of parameters and saving the results. sweep_config = fitting.SweepConfig( binary_sweep_parameters = binary_sweep_parameters, weighted_sweep_parameters = weighted_sweep_parameters, num_simulations = num_simulations, distance_matrix = [distance_matrix], )

In [ ]:

# Create our evaluation critereon - this tells you how close the generative model 
# is to your actual connectome in terms of a given set of topological or topographical 
# criteria - in this case, we use the KS Statistic to compare clustering, degree, and edge length

# binary evaluations
criteria = [ evaluation.ClusteringKS(), evaluation.DegreeKS(), evaluation.EdgeLengthKS(distance_matrix) ]
energy_equation = evaluation.MaxCriteria( criteria )
binary_evaluations = [energy_equation]

# weighted evaluations
weighted_evaluations = [ evaluation.WeightedNodeStrengthKS(normalise=True), evaluation.WeightedClusteringKS() ]

# Create our evaluation critereon - this tells you how close the generative model # is to your actual connectome in terms of a given set of topological or topographical # criteria - in this case, we use the KS Statistic to compare clustering, degree, and edge length # binary evaluations criteria = [ evaluation.ClusteringKS(), evaluation.DegreeKS(), evaluation.EdgeLengthKS(distance_matrix) ] energy_equation = evaluation.MaxCriteria( criteria ) binary_evaluations = [energy_equation] # weighted evaluations weighted_evaluations = [ evaluation.WeightedNodeStrengthKS(normalise=True), evaluation.WeightedClusteringKS() ]

Perform the sweep

In [ ]:

# so we can see how long the sweep takes, we can start a timer
# and stop it when the sweep is finished
start_time = time.perf_counter()

# Given all our parameter configurations, we can now run the sweep.
# This will take a minute, so be patient.
experiments = fitting.perform_sweep(
    sweep_config=sweep_config, 
    binary_evaluations=binary_evaluations, 
    real_binary_matrices=binary_consensus_network,
    real_weighted_matrices=weighted_consensus_network,
    weighted_evaluations=weighted_evaluations,
    save_model = False,
    save_run_history = False, 
    verbose=True # set this to True if you want to see a progress bar below
)

# end the timer
end_time = time.perf_counter()

# so we can see how long the sweep takes, we can start a timer # and stop it when the sweep is finished start_time = time.perf_counter() # Given all our parameter configurations, we can now run the sweep. # This will take a minute, so be patient. experiments = fitting.perform_sweep( sweep_config=sweep_config, binary_evaluations=binary_evaluations, real_binary_matrices=binary_consensus_network, real_weighted_matrices=weighted_consensus_network, weighted_evaluations=weighted_evaluations, save_model = False, save_run_history = False, verbose=True # set this to True if you want to see a progress bar below ) # end the timer end_time = time.perf_counter()

Configuration Iterations: 100%|██████████| 9/9 [01:17<00:00,  8.58s/it]

In [ ]:

# now we can see how long the sweep took

print(f"Sweep took {end_time - start_time:0.3f} seconds.")

total_simulations = num_simulations * len(alpha_values) * len(omega_values)

print(f"Total number of simulations: {total_simulations}")

print(f"Average time per simulation: {(end_time - start_time) / total_simulations:0.3f} seconds.")

# now we can see how long the sweep took print(f"Sweep took {end_time - start_time:0.3f} seconds.") total_simulations = num_simulations * len(alpha_values) * len(omega_values) print(f"Total number of simulations: {total_simulations}") print(f"Average time per simulation: {(end_time - start_time) / total_simulations:0.3f} seconds.")

Sweep took 77.251 seconds.
Total number of simulations: 900
Average time per simulation: 0.086 seconds.