Merge branch 'main' into mainfunction

Author: aajayi-21

diffpy/snmf/io.py

@@ -4,7 +4,7 @@
 from diffpy.utils.parsers.loaddata import loadData
 
 
-def initialize_variables(data_input, component_amount, data_type, sparsity=1, smoothness=1e18):
+def initialize_variables(data_input, number_of_components, data_type, sparsity=1, smoothness=1e18):
     """Determines the variables and initial values used in the SNMF algorithm.
 
     Parameters
@@ -13,7 +13,7 @@ def initialize_variables(data_input, component_amount, data_type, sparsity=1, sm
       The observed or simulated PDF or XRD data provided by the user. Has dimensions R x N where R is the signal length
       and N is the number of PDF/XRD signals.
 
-    component_amount: int
+    number_of_components: int
       The number of component signals the user would like to decompose 'data_input' into.
 
     data_type: str
@@ -40,27 +40,19 @@ def initialize_variables(data_input, component_amount, data_type, sparsity=1, sm
 
     """
     signal_length = data_input.shape[0]
-    moment_amount = data_input.shape[1]
+    number_of_signals = data_input.shape[1]
 
-    component_matrix_guess = np.random.rand(signal_length, component_amount)
-    weight_matrix_guess = np.random.rand(component_amount, moment_amount)
-    stretching_matrix_guess = np.ones(component_amount, moment_amount) + np.random.randn(component_amount,
-                                                                                         moment_amount) * 1e-3
+    diagonals = [np.ones(number_of_signals - 2), -2 * np.ones(number_of_signals - 2), np.ones(number_of_signals - 2)]
+    smoothness_term = .25 * scipy.sparse.diags(diagonals, [0, 1, 2], shape=(number_of_signals - 2, number_of_signals))
 
-    diagonals = [np.ones(moment_amount - 2), -2 * np.ones(moment_amount - 2), np.ones(moment_amount - 2)]
-    smoothness_term = .25 * scipy.sparse.diags(diagonals, [0, 1, 2], shape=(moment_amount - 2, moment_amount))
-
-    hessian_helper_matrix = scipy.sparse.block_diag([smoothness_term.T @ smoothness_term] * component_amount)
-    sequence = np.arange(moment_amount * component_amount).reshape(component_amount, moment_amount).T.flatten()
+    hessian_helper_matrix = scipy.sparse.block_diag([smoothness_term.T @ smoothness_term] * number_of_components)
+    sequence = np.arange(number_of_signals * number_of_components).reshape(number_of_components, number_of_signals).T.flatten()
     hessian_helper_matrix = hessian_helper_matrix[sequence, :][:, sequence]
 
     return {
         "signal_length": signal_length,
-        "moment_amount": moment_amount,
-        "component_matrix_guess": component_matrix_guess,
-        "weight_matrix_guess": weight_matrix_guess,
-        "stretching_matrix_guess": stretching_matrix_guess,
-        "component_amount": component_amount,
+        "number_of_signals": number_of_signals,
+        "number_of_components": number_of_components,
         "data_type": data_type,
         "smoothness": smoothness,
         "sparsity": sparsity,

diffpy/snmf/subroutines.py

@@ -25,7 +25,7 @@ def initialize_components(number_of_components, number_of_signals, grid_vector):
         raise ValueError(f"Number of components = {number_of_components}. Number_of_components must be >= 1.")
     components = list()
     for component in range(number_of_components):
-        component = ComponentSignal(grid_vector,number_of_signals,component)
+        component = ComponentSignal(grid_vector, number_of_signals, component)
         components.append(component)
     return tuple(components)
 
@@ -54,6 +54,36 @@ def lift_data(data_input, lift=1):
     return data_input + np.abs(np.min(data_input) * lift)
 
 
+def construct_stretching_matrix(components, number_of_components, number_of_signals):
+    """Constructs the stretching factor matrix
+
+    Parameters
+    ----------
+    components: tuple of ComponentSignal objects
+       The tuple containing the component signals in ComponentSignal objects.
+    number_of_signals: int
+      The number of signals in the data provided by the user.
+
+    Returns
+    -------
+    2d array
+      The matrix containing the stretching factors for the component signals for each of the signals in the raw data.
+      Has dimensions `component_signal` x `number_of_signals`
+
+    """
+    if (len(components)) == 0:
+        raise ValueError(f"Number of components = {number_of_components}. Number_of_components must be >= 1.")
+    number_of_components = len(components)
+    
+    if number_of_signals <= 0:
+        raise ValueError(f"Number of signals = {number_of_signals}. Number_of_signals must be >= 1.")
+        
+    stretching_factor_matrix = np.zeros((number_of_components, number_of_signals))
+    for i, component in enumerate(components):
+        stretching_factor_matrix[i] = component.stretching_factors
+    return stretching_factor_matrix
+
+
 def initialize_arrays(number_of_components, number_of_moments, signal_length):
     """Generates the initial guesses for the weight, stretching, and component matrices
 

diffpy/snmf/tests/test_subroutines.py

@@ -1,7 +1,8 @@
 import pytest
 import numpy as np
+from diffpy.snmf.containers import ComponentSignal
 from diffpy.snmf.subroutines import objective_function, get_stretched_component, reconstruct_data, get_residual_matrix, \
-    update_weights_matrix, initialize_arrays, lift_data, initialize_components
+    update_weights_matrix, initialize_arrays, lift_data, initialize_components, construct_stretching_matrix
 
 to = [
     ([[[1, 2], [3, 4]], [[5, 6], [7, 8]], 1e11, [[1, 2], [3, 4]], [[1, 2], [3, 4]], 1], 2.574e14),
@@ -161,3 +162,27 @@ def test_initialize_components(tcc):
     assert len(actual) == tcc[0]
     assert len(actual[0].weights) == tcc[1]
     assert (actual[0].grid == np.array(tcc[2])).all()
+
+tcso =[([ComponentSignal([0,.5,1,1.5],20,0)],1,20),
+       ([ComponentSignal([0,.5,1,1.5],20,0)],4,20),
+       # ([ComponentSignal([0,.5,1,1.5],20,0)],0,20), # Raises an exception
+       # ([ComponentSignal([0,.5,1,1.5],20,0)],-2,20), # Raises an exception
+       # ([ComponentSignal([0,.5,1,1.5],20,0)],1,0), # Raises an Exception
+       # ([ComponentSignal([0,.5,1,1.5],20,0)],1,-3), # Raises an exception
+       ([ComponentSignal([0,.5,1,1.5],20,0),ComponentSignal([0,.5,1,1.5],20,1)],2,20),
+       ([ComponentSignal([0,.5,1,1.5],20,0),ComponentSignal([0,.5,1,21.5],20,1)],2,20),
+       ([ComponentSignal([0,1,1.5],20,0),ComponentSignal([0,.5,1,21.5],20,1)],2,20),
+       # ([ComponentSignal([0,.5,1,1.5],20,0),ComponentSignal([0,.5,1,1.5],20,1)],1,-3), # Negative signal length. Raises an exception
+       #([],1,20), # Empty components. Raises an Exception
+       #([],-1,20), # Empty components with negative number of components. Raises an exception
+       #([],0,20), # Empty components with zero number of components. Raises an exception
+       #([],1,0), # Empty components with zero signal length. Raises an exception.
+       #([],-1,-2), # Empty components with negative number of components and signal length Raises an exception.
+
+]
+@pytest.mark.parametrize('tcso',tcso)
+def test_construct_stretching_matrix(tcso):
+    actual = construct_stretching_matrix(tcso[0],tcso[1],tcso[2])
+    for component in tcso[0]:
+        np.testing.assert_allclose(actual[component.id,:], component.stretching_factors)
+        #assert actual[component.id, :] == component.stretching_factors

setup.py

@@ -1,6 +1,6 @@
 #!/usr/bin/env python
 
-# Installation script for diffpy.pdfmorph
+# Installation script for diffpy.snmf
 
 """diffpy.snmf - a package implementing the stretched NMF algorithm.
 
@@ -28,7 +28,7 @@
     entry_points={
         # define console_scripts here, see setuptools docs for details.
         'console_scripts': [
-            'pdfmorph = diffpy.pdfmorph.pdfmorphapp:main',
+            'snmf = diffpy.snmf.stretchednmfapp:main',
         ],
     },
     test_suite='tests',