From 57a4dfa94edb4272a64363094faa81649b5b3991 Mon Sep 17 00:00:00 2001
From: Blair Lyons <blair208@gmail.com>
Date: Tue, 3 May 2022 14:27:21 -0700
Subject: [PATCH] add metrics from ActinAnalyzer

---
 .../processes/readdy_actin_process.py         | 77 ++++++++++++++++++-
 vivarium_models/util.py                       | 24 ++++++
 2 files changed, 100 insertions(+), 1 deletion(-)

diff --git a/vivarium_models/processes/readdy_actin_process.py b/vivarium_models/processes/readdy_actin_process.py
index b69e4e4..e6a26df 100644
--- a/vivarium_models/processes/readdy_actin_process.py
+++ b/vivarium_models/processes/readdy_actin_process.py
@@ -6,9 +6,15 @@
 from vivarium.core.control import run_library_cli
 
 from tqdm import tqdm
-from simularium_models_util.actin import ActinSimulation, ActinUtil, ActinTestData
+from simularium_models_util.actin import (
+    ActinSimulation,
+    ActinUtil,
+    ActinTestData,
+    ActinAnalyzer,
+)
 from simularium_models_util import ReaddyUtil
 from vivarium_models.util import create_monomer_update
+from vivarium_models.util import create_monomer_update, format_monomer_results
 from vivarium_models.library.scan import Scan
 
 NAME = "ReaDDy_actin"
@@ -318,10 +324,79 @@ def filament_lengths(results):
         length = np.linalg.norm(difference)
         return length
 
+    def percent_filamentous_actin(results):
+        monomer_data = [list(results.values())[-1]["monomers"]]
+        return (
+            100.0
+            * ActinAnalyzer.analyze_ratio_of_filamentous_to_total_actin(monomer_data)[0]
+        )
+
+    def percent_bound_arp23(results):
+        monomer_data = [list(results.values())[-1]["monomers"]]
+        return (
+            100.0 * ActinAnalyzer.analyze_ratio_of_bound_to_total_arp23(monomer_data)[0]
+        )
+
+    def percent_daughter_actin(results):
+        monomer_data = [list(results.values())[-1]["monomers"]]
+        return (
+            100.0
+            * ActinAnalyzer.analyze_ratio_of_daughter_to_total_actin(monomer_data)[0]
+        )
+
+    def mother_filament_lengths(results):
+        monomer_data = [list(results.values())[-1]["monomers"]]
+        return ActinAnalyzer.analyze_average_over_time(
+            ActinAnalyzer.analyze_mother_filament_lengths(monomer_data)
+        )[0]
+
+    def daughter_filament_lengths(results):
+        monomer_data = [list(results.values())[-1]["monomers"]]
+        return ActinAnalyzer.analyze_daughter_filament_lengths(monomer_data)[0]
+
+    def branch_angles(results):
+        box_size = list(results.values())[-1]["monomers"]["box_size"]
+        monomer_data = [list(results.values())[-1]["monomers"]]
+        periodic_boundary = False  # TODO get from parameters
+        return ActinAnalyzer.analyze_branch_angles(
+            monomer_data, box_size, periodic_boundary
+        )[0]
+
+    def short_helix_pitches(results):
+        monomer_data = format_monomer_results(results)
+        periodic_boundary = False  # TODO get from parameters
+        return ActinAnalyzer.analyze_short_helix_pitches(
+            monomer_data, monomer_data[0]["box_size"], periodic_boundary
+        )
+
+    def long_helix_pitches(results):
+        monomer_data = format_monomer_results(results)
+        periodic_boundary = False  # TODO get from parameters
+        return ActinAnalyzer.analyze_long_helix_pitches(
+            monomer_data, monomer_data[0]["box_size"], periodic_boundary
+        )
+
+    def filament_straightness(results):
+        monomer_data = format_monomer_results(results)
+        periodic_boundary = False  # TODO get from parameters
+        return ActinAnalyzer.analyze_filament_straightness(
+            monomer_data, monomer_data[0]["box_size"], periodic_boundary
+        )
+
     metrics = {
         'count_monomers': count_monomers,
         'count_monomer_types': count_monomer_types,
         'filament_lengths': filament_lengths}
+        "percent_filamentous_actin": percent_filamentous_actin,
+        "percent_bound_arp23": percent_bound_arp23,
+        "percent_daughter_actin": percent_daughter_actin,
+        "mother_filament_lengths": mother_filament_lengths,
+        "daughter_filament_lengths": daughter_filament_lengths,
+        "branch_angles": branch_angles,
+        "short_helix_pitches": short_helix_pitches,
+        "long_helix_pitches": long_helix_pitches,
+        "filament_straightness": filament_straightness,
+    }
 
     scan = Scan(
         parameters,
diff --git a/vivarium_models/util.py b/vivarium_models/util.py
index 1c5341e..9be67e5 100644
--- a/vivarium_models/util.py
+++ b/vivarium_models/util.py
@@ -1,3 +1,6 @@
+import numpy as np
+
+
 def agents_update(existing, projected):
     update = {"_add": [], "_delete": []}
 
@@ -26,3 +29,24 @@ def create_monomer_update(previous_monomers, new_monomers):
     return {
         "monomers": {"topologies": topologies_update, "particles": particles_update}
     }
+
+
+def format_monomer_results(results):
+    """
+    Workaround since numpy arrays are not preserved in Vivarium?
+    """
+    monomer_data = list(results.values())[-1]["monomers"]
+    formatted_results = {
+        "box_center": monomer_data["box_center"],
+        "box_size": monomer_data["box_size"],
+        "topologies": monomer_data["topologies"],
+        "particles": {},
+    }
+    for particle_id in monomer_data["particles"]:
+        particle = monomer_data["particles"][particle_id]
+        formatted_results["particles"][particle_id] = {
+            "type_name": particle["type_name"],
+            "position": np.array(particle["position"]),
+            "neighbor_ids": particle["neighbor_ids"],
+        }
+    return [formatted_results]