Merge pull request #17 from s-ccs/fix-sparseInit

fix docs simulation, had overlap where it shouldnt have had any, also sometimes sparse matrix was a bit too short

Author: behinger

docs/src/reference/type1.jl

@@ -41,20 +41,20 @@ function run_fun(r)
         MersenneTwister(r),
         design,
         signal,
-        UniformOnset(; offset=5, width=4),
+        NoOnset(),
         RedNoise(noiselevel=1);
         return_epoched=true,
     )
     data = reshape(data, size(data, 1), :)
-    data = data[:, evts.condition.=="small"] .- data[:, evts.condition.=="large"]
+    data_diff = data[:, evts.condition.=="small"] .- data[:, evts.condition.=="large"]
 
     return data,
-    clusterdepth(data'; τ=quantile(TDist(n_subjects - 1), 0.975), nperm=1000, show_warnings=false)
+    clusterdepth(data_diff; τ=quantile(TDist(n_subjects - 1), 0.975), nperm=4000, show_warnings=false)
 end;
 
 # ## Understanding the simulation
 # let's have a look at the actual data by running it once, plotting condition wise trials, the ERP and histograms of uncorrected and corrected p-values
-data, pval = run_fun(5)
+data, pval = run_fun(1)
 conditionSmall = data[:, 1:2:end]
 conditionLarge = data[:, 2:2:end]
 pval_uncorrected =
@@ -96,7 +96,7 @@ reps = 500
 res = fill(NaN, reps, 2)
 Threads.@threads for r = 1:reps
     data, pvals = run_fun(r)
-    res[r, 1] = mean(pvals .<= 0.05 / 2)
+    res[r, 1] = mean(pvals .<= 0.05)
     res[r, 2] =
         mean(abs.(ClusterDepth.studentt(data)) .>= quantile(TDist(n_subjects - 1), 0.975))
 end;

src/cluster.jl

@@ -105,9 +105,11 @@ function perm_clusterdepths_both(
     end
     #@debug size(d0)
     #@debug size(data_perm)
+    #data_perm .-= mean(data_perm, dims=length(size(data_perm)))
     for i = 1:nₚ
         # permute	
         d_perm = permfun(rng, data_perm)
+
         if isnothing(statfun!)
             d0 = statfun(d_perm)
         else
@@ -136,8 +138,8 @@ function perm_clusterdepths_both(
         end
     end
 
-    Jₖ_head = sparse(rows_h, cols_h, vals_h)#SparseMatrixCSC(nₚ,maximum(rows_h), cols_h,rows_h,vals_h)
-    Jₖ_tail = sparse(rows_t, cols_t, vals_t)#SparseMatrixCSC(nₚ,maximum(rows_t), cols_t,rows_t,vals_t)
+    Jₖ_head = sparse(rows_h, cols_h, vals_h, maximum(rows_h), nₚ)#SparseMatrixCSC(nₚ,maximum(rows_h), cols_h,rows_h,vals_h)
+    Jₖ_tail = sparse(rows_t, cols_t, vals_t, maximum(rows_t), nₚ)#SparseMatrixCSC(nₚ,maximum(rows_t), cols_t,rows_t,vals_t)
     return ClusterDepthMatrix((Jₖ_head)), ClusterDepthMatrix((Jₖ_tail))
 end
 
@@ -264,3 +266,4 @@ function cluster(data::BitVector)
     end
     return start, len
 end
+show_warnings

test/cluster.jl

@@ -1,3 +1,5 @@
+
+
 @testset "cluster" begin
     s, l = ClusterDepth.cluster(
         [4.0, 0.0, 10.0, 0.0, 3.0, 4.0, 0, 4.0, 4.0, 0.0, 0.0, 5.0] .> 0.9,
@@ -64,4 +66,12 @@ end
     data[23, :] .-= 3
     @test_warn x -> occursin("Your data shows a cluster", x) ClusterDepth.clusterdepth(data; τ=0.4, nperm=5)
 
+end
+
+@testset "sparse matrix size" begin
+    data = rand(StableRNG(42), 5, 10) .- 0.5
+    data[2, :] .+= 1
+    data[3, :] .+= 3
+    cdmTuple = ClusterDepth.perm_clusterdepths_both(StableRNG(1), data, (ClusterDepth.sign_permute!), 1.5; nₚ=1000, (statfun!)=(ClusterDepth.studentt!), sidefun=abs)
+    @test size(cdmTuple[1].J, 2) == 1000 # sparse matrix has to have the size of n-perms
 end