threads limitatiton and "Error in if (sol_epgap > opt$epgap) { : the condition has length > 1"

[9421g]

Hello, This error message is from jabba with R 4.3.3 and gurobi. Is this a bug of jabba?

Another problem is Gurobi is using all cores instead of my setting jba --cores 3. I've tried using echo "Threads 1" > gurobi.env, but It seems not working.

https://support.gurobi.com/hc/en-us/articles/14126085438481-How-do-I-set-parameters-in-Gurobi

It's ends up that gGnome::run_gurobi did not set threads limitation:

    ## params
    params = list()
    if (!is.null(control$epgap)) {
        params$MIPGap = control$epgap
    }
    if (!is.null(control$tilim)) {
        params$TimeLimit = control$tilim
    }
    if (!is.null(control$trace)) {
        params$LogToConsole = ifelse(control$trace > 0, 1, 0)
    }

    ## TODO: set up env list for running on compute cluster

    ## run gurobi
    sol = gurobi::gurobi(model = model, params = params)

And in jbaMIP:

sol = gurobi::gurobi(model, params = c(list(TimeLimit=tilim), list(...)));

I turned this packages into a script and solved these anyway.

detected no tier 1 junctions
No duplicate Tier 1 junctions detected
Grabbing available memory...
Currently used: 1132.1 Mb
Allowed: 16000 Mb
Treemem: 13867.9 Mb
creating copy of input gGraph
Checking inputs
adding l0 penalty indicator
adding delta constraints for LP
Unique cids (A): 377724
Unique cids (b): 377724
Number of variables: 374440
bvec length: 377724
Amat nrow: 377724
Starting optimization with gurobi!
Gurobi Optimizer version 10.0.3 build v10.0.3rc0 (linux64)

CPU model: Intel(R) Xeon(R) Gold 5218R CPU @ 2.10GHz, instruction set [SSE2|AVX|AVX2|AVX512]
Thread count: 40 physical cores, 40 logical processors, using up to 32 threads

Optimize a model with 377724 rows, 374444 columns and 808884 nonzeros
Model fingerprint: 0x9a945994
Variable types: 180492 continuous, 193952 integer (66896 binary)
Coefficient statistics:
  Matrix range     [1e-01, 1e+03]
  Objective range  [6e-01, 4e+05]
  Bounds range     [1e+03, 1e+03]
  RHS range        [8e-02, 2e+01]
Found heuristic solution: objective 2.587947e+07
Presolve removed 409025 rows and 444930 columns (presolve time = 9s) ...
Presolve removed 368988 rows and 363824 columns
Presolve time: 9.41s
Presolved: 8736 rows, 10620 columns, 28087 nonzeros
Found heuristic solution: objective 2873465.6700
Variable types: 3625 continuous, 6995 integer (3498 binary)

Root simplex log...

Iteration    Objective       Primal Inf.    Dual Inf.      Time
       0    7.2680484e+05   4.470832e+02   0.000000e+00     10s
    4814    7.2716900e+05   0.000000e+00   0.000000e+00     10s

Root relaxation: objective 7.271690e+05, 4814 iterations, 0.03 seconds (0.01 work units)

    Nodes    |    Current Node    |     Objective Bounds      |     Work
 Expl Unexpl |  Obj  Depth IntInf | Incumbent    BestBd   Gap | It/Node Time

     0     0 727169.004    0 4340 2873465.67 727169.004  74.7%     -   10s
H    0     0                    1219078.8603 727169.004  40.4%     -   10s
H    0     0                    1082296.1985 727169.004  32.8%     -   10s
H    0     0                    1079457.7230 727169.004  32.6%     -   10s
H    0     0                    887061.11027 727169.004  18.0%     -   10s
H    0     0                    883049.78607 727169.004  17.7%     -   10s
     0     0 760787.893    0 2877 883049.786 760787.893  13.8%     -   10s
H    0     0                    870402.54721 760787.893  12.6%     -   10s
H    0     0                    865956.97232 760787.893  12.1%     -   10s
H    0     0                    850238.83050 760787.893  10.5%     -   10s
H    0     0                    821108.89250 785052.588  4.39%     -   10s
     0     0 785056.039    0 2052 821108.892 785056.039  4.39%     -   10s
     0     0 785069.364    0 2059 821108.892 785069.364  4.39%     -   10s
H    0     0                    816524.74103 813927.750  0.32%     -   10s
     0     0 814026.919    0  105 816524.741 814026.919  0.31%     -   10s
H    0     0                    816145.71175 814026.919  0.26%     -   10s
H    0     0                    815006.17743 814161.371  0.10%     -   10s
     0     0 814161.371    0   47 815006.177 814161.371  0.10%     -   10s
     0     0 814266.988    0   28 815006.177 814266.988  0.09%     -   10s
H    0     0                    814706.00917 814266.988  0.05%     -   11s
H    0     0                    814616.95537 814266.988  0.04%     -   11s
H    0     0                    814300.38399 814266.988  0.00%     -   11s
     0     0 814299.273    0    9 814300.384 814299.273  0.00%     -   11s

Cutting planes:
  Gomory: 197
  Implied bound: 99
  MIR: 2993
  StrongCG: 1
  Flow cover: 13
  RLT: 21
  Relax-and-lift: 81

Explored 1 nodes (8408 simplex iterations) in 11.20 seconds (5.87 work units)
Thread count was 32 (of 40 available processors)

Solution count 10: 814300 814617 814706 ... 870403

Optimal solution found (tolerance 1.00e-05)
Best objective 8.143003839902e+05, best bound 8.142992729000e+05, gap 0.0001%
CPLEX epgap 1.3644720798363e-06 with solution status 
JaBbA 2024-11-14 11:09:42.829018: Recording convergence status of subgraphs
JaBbA 2024-11-14 11:09:50.792997: simplifying segments in JaBbA graph but keeping all edges (including copy 0), dumping to jabba.rds
JaBbA 2024-11-14 11:10:00.536339: Checking for hets
JaBbA 2024-11-14 11:10:00.800193: Skipping loose end annotation
JaBbA 2024-11-14 11:10:00.800695: Saving results
JaBbA 2024-11-14 11:10:13.249471: Generating figures
JaBbA 2024-11-14 11:10:42.222631: Done .. job output in: /data/person/cells/taokj/chenja/chromothripsis/3.gGnome/NANT01002-T0-3A
Error in if (sol_epgap > opt$epgap) { : the condition has length > 1
Calls: suppressPackageStartupMessages -> withCallingHandlers -> JaBbA
In addition: There were 26 warnings (use warnings() to see them)
Execution halted

[shihabdider]

So this error arises from a final guard which checks if the solution epgap meets the threshold. I would double check the opt$epgap being passed, because this error indicates that it's passed as multiple values instead of a single number.

[9421g]

Hi @shihabdider, thank you for reply. Another bug occured. I run jabba on 71 samples, 57 successfully output jabba.simple.gg.rds, 14 failed with:

Error in orimap[sapply(seq_along(ori), function(i) ori[[i]][iid[i]])] : 
  invalid subscript type 'list'
Calls: suppressPackageStartupMessages -> withCallingHandlers -> JaBbA -> read.junctions
In addition: Warning messages:
1: JaBbA 2024-11-16 04:13:15.006408: Vcf not in proper format.  Will treat rearrangements as if in BND format 
2: JaBbA 2024-11-16 04:13:15.012265: 602 rows of vcf do not have svtype BND, treat them as non-BND! 

My junction input file is Delly raw output. I simply solved this by using

rtracklayer::export(
    StructuralVariantAnnotation::breakpointgr2pairs(
        StructuralVariantAnnotation::partner(
            StructuralVariantAnnotation::breakpointRanges(
                VariantAnnotation::readVcf(
                        commandArgs(trailingOnly = TRUE)[1]
                ),
                inferMissingBreakends = TRUE
            )
        )
    ),
    con = commandArgs(trailingOnly = TRUE)[2]
)

to convert vcf to bedpe.

[9421g]

And this is another bug:

JaBbA 2024-11-16 13:25:14.650338: Number of gaps with nonzero width: 6732
JaBbA 2024-11-16 13:25:14.650628: Number of segments before gap filtering: 27165
JaBbA 2024-11-16 13:25:15.206296: 20433 segments produced after gap filtering
JaBbA 2024-11-16 13:25:17.242545: Brand new function for reciprocal junctions calling.
..
JaBbA 2024-11-16 13:25:18.447472: Excluding 1 aberrant junctions whose both breakpoints are in NA coverage regions
JaBbA 2024-11-16 13:25:18.448: Cancel nudge for 18 aberrant junctions where one of the 2 breakpoint is in NA coverage regions
JaBbA 2024-11-16 13:25:18.448293: In sum, we are forcing 0 junctions, excluding 1 junctions, and nudging 14 junctions
JaBbA 2024-11-16 13:25:19.293571: Excluding aberrant junctions: 509
JaBbA 2024-11-16 13:25:19.377128: Adjusting the kag (naive solution) as mipstart (initial solution).
Error in `:=`(cn, cnmle) : 
  Check that is.data.table(DT) == TRUE. Otherwise, :=, `:=`(...) and let(...) are defined for use in j, once only and in particular ways. Note that namespace-qualification like data.table::`:=`(...) is not supported. See help(":=").
Calls: suppressPackageStartupMessages ... [.data.frame -> := -> stopf -> raise_condition -> signal
In addition: Warning messages:
1: JaBbA 2024-11-16 13:24:39.904294: Tier field tier missing: giving every junction the same tier, i.e. all have the potential to be incorporated 
2: JaBbA 2024-11-16 13:25:15.248105: Skipping over karyograph creation because file already exists and overwrite = FALSE 
3: `clusters()` was deprecated in igraph 2.0.0.
ℹ Please use `components()` instead. 
4: `graph.adjacency()` was deprecated in igraph 2.0.0.
ℹ Please use `graph_from_adjacency_matrix()` instead. 
5: In `[.data.table`(junc.dt, , `:=`(from = NULL, to = NULL)) :
  Tried to assign NULL to column 'from', but this column does not exist to remove
6: In `[.data.table`(junc.dt, , `:=`(from = NULL, to = NULL)) :
  Tried to assign NULL to column 'to', but this column does not exist to remove
Execution halted