Before christmas

clustering/class/dissimilarity_matrix.py

@@ -4,6 +4,7 @@
 import numpy.ma as ma
 
 from numba import jit
+from scipy.sparse import csr_matrix
 from scipy.spatial.distance import pdist, squareform
 from sklearn.metrics import pairwise_distances
 from sklearn.preprocessing import MinMaxScaler
@@ -157,7 +158,7 @@ def diagonal_dist(self, u, v):
         d2 = abs(v[0] - v[1])
 
         dist = abs(d1 - d2)/self.n
-        return dist
+        return np.log10(dist)
 
     def manhattan(self, u, v):
         """Calculates the scaled manahattan distance between u and v.
@@ -223,7 +224,7 @@ def sklearn_dist(self, metric, scaler=None, **kwargs):
 
         # calculate the pairwise distances between data point and convert it
         # to square distance matrix
-        distmat = pairwise_distances(self.X, metric=metric, force_all_finite=False, 
+        distmat = pairwise_distances(csr_matric(self.X), metric=metric, force_all_finite=False, 
                                      n_jobs=-1, **kwargs)
 
         return distmat

clustering/class/simulate_data.py

@@ -2,6 +2,37 @@
 
 from data_preparation import DataPreparation
 
+def _powerlaw(distance, C, alpha):
+    result = np.where(distance == 0, C, C * distance**(-alpha))
+    return result
+
+def _poisson(distance):
+    result = np.zeros(distance)
+
+    for i, lmbda in enumerate(distance):
+        result[i] = np.random.Generator.poisson(lmbda)
+    return result
+
+def _sim_mat(n, medianIF, powerlaw_alpha):
+    cell1 = np.zeros((n, n))
+    cell2 = np.zeros((n, n))
+
+    ij = np.ones((n, n))*np.linspace(0, n-1, n)
+    i = ij.T[np.triu_indices(n)]
+    j = ij[np.triu_indices(n)]
+
+    distance = j - i + 1
+
+    #Bmean = _powerlaw(distance, medianIF, powerlaw_alpha)
+    Bmean = _poisson(distance)
+
+    cell1[np.triu_indices(n)] = Bmean
+    cell1 = cell1 + cell1.T
+    cell1[np.diag_indices(n)] /= 2
+
+    return cell1
+
+
 class SimulateData():
     def __init__(self, resolution, region):
         """
@@ -43,10 +74,19 @@ def initial_matrix(self, resolution, region):
         wt_001 = DataPreparation(path_wt_001, resolution, region)
         wt_002 = DataPreparation(path_wt_002, resolution, region)
 
-        matrix = wt_002.subtract(wt_001)
+        matrix = wt_001.matrix - wt_002.subtract(wt_001)
 
         return wt_001.clr, matrix
 
+    def addTAD(self, tad_mat, ul):
+        test = np.zeros_like(self.matrix)
+        test[ul:ul + tad_mat.shape[0], ul:ul + tad_mat.shape[1]] = tad_mat
+
+        return self.matrix + test
+
+
+
+
 
 
 

clustering/class/untitled.py

@@ -0,0 +1,11 @@
+import numpy as np
+
+class ContactMatrix:
+    def __init__(self):
+        pass
+
+    def transpose_coords(self, row, col):
+        if row > col:
+            return (row - col, row)
+        else:
+            return (col - row, col)

clustering/notebooks/simulate_hic.py

@@ -0,0 +1,24 @@
+import numpy as np
+
+def tads(k, distance, alpha):
+    interaction = k*distance**alpha
+    return interaction
+
+def _sim_mat(n, k, alpha):
+    cell1 = np.zeros((n, n))
+
+    ij = np.ones((n, n))*np.linspace(0, n-1, n)
+    i = ij.T[np.triu_indices(n)]
+    j = ij[np.triu_indices(n)]
+
+    distance = j - i + 1
+
+    Bmean = tads(k, distance, alpha)
+
+    cell1[np.triu_indices(n)] = Bmean
+    cell1 = cell1 + cell1.T
+    cell1[np.diag_indices(n)] /= 2
+
+    return cell1
+
+class SimulateHiC:

examples/helloworld.rs

@@ -0,0 +1,3 @@
+fn main() {
+    println!("Hello World!");
+}

examples/test.rs

@@ -0,0 +1,8 @@
+extern crate rand;
+
+fn main() {
+let dissim = ndarray::arr2(&[[0,1,2,3],[1,0,4,5],[2,4,0,6],[3,5,6,0]]);
+let mut meds = kmedoids::random_initialization(4, 2, &mut rand::thread_rng());
+let (loss, assingment, n_iter, n_swap) = kmedoids::fasterpam(&dissim, &mut meds, 100);
+println!("Loss is: {}", loss);
+}

examples/testing/Cargo.toml

@@ -0,0 +1,9 @@
+[package]
+name = "testing"
+version = "0.1.0"
+edition = "2021"
+
+# See more keys and their definitions at https://doc.rust-lang.org/cargo/reference/manifest.html
+
+[dependencies]
+rand = "0.3.14"

examples/testing/src/main.rs

@@ -0,0 +1,3 @@
+fn main() {
+    println!("Hello, world!");
+}

examples/testing/test_rand.rs

@@ -0,0 +1,14 @@
+use rand::Rng;
+
+fn main() {
+    let mut rng = rand::thread_rng();
+
+    let n1: u8 = rng.gen();
+    let n2: u16 = rng.gen();
+    println!("Random u8: {}", n1);
+    println!("Random u16: {}", n2);
+    println!("Random u32: {}", rng.gen::<u32>());
+    println!("Random i32: {}", rng.gen::<i32>());
+    println!("Random float: {}", rng.gen::<f64>());
+}
+