AlphaGenes / AlphaPeel

AlphaPeel: calling, phasing, and imputing genotype and sequence data in pedigrees
MIT License
2 stars 11 forks source link

Sort out numba issue: 'Don't know how to allocate array with layout 'A'.' #7

Closed gregorgorjanc closed 1 year ago

gregorgorjanc commented 1 year ago

I am experiencing issue described by https://stackoverflow.com/questions/74681366/alphapeel-numba-fails-to-allocate-array-with-layout-a, which seems to be the same as https://stackoverflow.com/questions/64044846/python-slicing-2d-numpy-array-to-produce-order-c-arrays-while-using-numba

I tried using childValues = posterior[child,:,:] * penetrance[child,:,:] or childValues = np.ascontiguousarray(posterior[child,:,:] * penetrance[child,:,:]), but both give the same error:(

Here is reproducible example code

git clone https://github.com/AlphaGenes/AlphaPeel.git
cd AlphaPeel
git submodule init
git submodule update

mamba create -n AlphaPeel
mamba activate AlphaPeel
./build_pipeline.sh
pip install --force-reinstall dist/AlphaPeel*.whl

cd example
AlphaPeel -genotypes data/genotypes.txt -pedigree data/pedigree.txt -out outputs/multilocus -nCycles 1 -runType multi -maxthreads 1  

which returns

Reading in AlphaGenes format: data/genotypes.txt
Cycle  0
Peeling Down, Generation 0
Peeling Down, Generation 1
Traceback (most recent call last):
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/errors.py", line 823, in new_error_context
    yield
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/lowering.py", line 265, in lower_block
    self.lower_inst(inst)
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/lowering.py", line 439, in lower_inst
    val = self.lower_assign(ty, inst)
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/lowering.py", line 626, in lower_assign
    return self.lower_expr(ty, value)
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/lowering.py", line 1368, in lower_expr
    res = self.context.special_ops[expr.op](self, expr)
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/np/ufunc/array_exprs.py", line 405, in _lower_array_expr
    return npyimpl.numpy_ufunc_kernel(
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/np/npyimpl.py", line 360, in numpy_ufunc_kernel
    output = _build_array(context, builder, ret_ty, sig.args, arguments)
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/np/npyimpl.py", line 286, in _build_array
    array_val = arrayobj._empty_nd_impl(context, builder, real_array_ty,
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/np/arrayobj.py", line 3923, in _empty_nd_impl
    raise NotImplementedError(
NotImplementedError: Don't know how to allocate array with layout 'A'.

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/bin/AlphaPeel", line 8, in <module>
    sys.exit(main())
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/tinypeel/tinypeel.py", line 255, in main
    runPeelingCycles(pedigree, peelingInfo, args, singleLocusMode = singleLocusMode)
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/tinypeel/tinypeel.py", line 24, in runPeelingCycles
    peelingCycle(pedigree, peelingInfo, args = args, singleLocusMode = singleLocusMode)
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/tinypeel/tinypeel.py", line 47, in peelingCycle
    Peeling.peel(family, Peeling.PEEL_DOWN, peelingInfo, singleLocusMode)
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/dispatcher.py", line 487, in _compile_for_args
    raise e
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/dispatcher.py", line 420, in _compile_for_args
    return_val = self.compile(tuple(argtypes))
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/dispatcher.py", line 965, in compile
    cres = self._compiler.compile(args, return_type)
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/dispatcher.py", line 125, in compile
    status, retval = self._compile_cached(args, return_type)
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/dispatcher.py", line 139, in _compile_cached
    retval = self._compile_core(args, return_type)
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/dispatcher.py", line 152, in _compile_core
    cres = compiler.compile_extra(self.targetdescr.typing_context,
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/compiler.py", line 716, in compile_extra
    return pipeline.compile_extra(func)
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/compiler.py", line 452, in compile_extra
    return self._compile_bytecode()
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/compiler.py", line 520, in _compile_bytecode
    return self._compile_core()
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/compiler.py", line 499, in _compile_core
    raise e
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/compiler.py", line 486, in _compile_core
    pm.run(self.state)
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/compiler_machinery.py", line 368, in run
    raise patched_exception
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/compiler_machinery.py", line 356, in run
    self._runPass(idx, pass_inst, state)
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/compiler_lock.py", line 35, in _acquire_compile_lock
    return func(*args, **kwargs)
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/compiler_machinery.py", line 311, in _runPass
    mutated |= check(pss.run_pass, internal_state)
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/compiler_machinery.py", line 273, in check
    mangled = func(compiler_state)
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/typed_passes.py", line 394, in run_pass
    lower.lower()
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/lowering.py", line 168, in lower
    self.lower_normal_function(self.fndesc)
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/lowering.py", line 222, in lower_normal_function
    entry_block_tail = self.lower_function_body()
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/lowering.py", line 251, in lower_function_body
    self.lower_block(block)
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/lowering.py", line 263, in lower_block
    with new_error_context('lowering "{inst}" at {loc}', inst=inst,
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/contextlib.py", line 153, in __exit__
    self.gen.throw(typ, value, traceback)
  File "/usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/numba/core/errors.py", line 837, in new_error_context
    raise newerr.with_traceback(tb)
numba.core.errors.LoweringError: Failed in nopython mode pipeline (step: native lowering)
Don't know how to allocate array with layout 'A'.

File "../../../../../../../../usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/tinypeel/Peeling/Peeling.py", line 86:
def peel(family, operation, peelingInfo, singleLocusMode) :
    <source elided>
        # We are estimating the parent's genotypes so the anterior term is ignored to avoid double counting.
        childValues = posterior[child,:,:] * penetrance[child,:,:]
        ^

During: lowering "childValues = arrayexpr(expr=(<built-in function mul>, [Var($528binary_subscr.17, Peeling.py:86), Var($548binary_subscr.29, Peeling.py:86)]), ty=array(float32, 2d, A))" at /usr/local/Caskroom/mambaforge/base/envs/AlphaPeel2/lib/python3.10/site-packages/tinypeel/Peeling/Peeling.py (86)
gregorgorjanc commented 1 year ago

Changing

@jit(nopython=True, nogil=True, locals={'e':float32, 'e4':float32, 'e16':float32, 'e1e':float32, 'childValues':float32[:,:]})

to

@jit(nopython=True, nogil=True, locals={'e': float32, 'e4':float32, 'e16':float32, 'e1e':float32})

makes the code run, and considerably faster

gregorgorjanc commented 1 year ago

Numba docs say The locals dictionary may be used to force the [Types and signatures](https://numba.pydata.org/numba-doc/latest/reference/types.html#numba-types) of particular local variables, for example if you want to force the use of single precision floats at some point. In general, we recommend you let Numba’s compiler infer the types of local variables by itself. https://numba.pydata.org/numba-doc/latest/reference/jit-compilation.html

gregorgorjanc commented 1 year ago

I think that indeed we don't need to state that childValues is float32[:,:] in @jit() in Peeling.py since it's created from

...
penetrance = peelingInfo.penetrance
posterior = peelingInfo.posterior
...
childValues = posterior[child,:,:] * penetrance[child,:,:]

where peelingInfo is created in class jit_peelingInformation(object): (in PeelingInfo.py) with these components defined as

...
self.posterior = np.full((self.nInd, 4, self.nLoci), baseValue, dtype = np.float32)
self.penetrance = np.full((self.nInd, 4, self.nLoci), baseValue, dtype = np.float32)
...

Also, other similar variables aren't specified in @jit() in Peeling.py, say:

probSire = anterior[sire,:,:]*penetrance[sire,:,:] * peelingInfo.posteriorSire_minusFam[fam,:,:]
gregorgorjanc commented 1 year ago

After implementing these changes I get

% Rscript checkAccuracy.r 
[1] " "
[1] "Assessing peeling file: multilocus.dosages"
[1] "Comparing accuracies:  0.98"
[1] " "
[1] "Assessing peeling file: hybrid.dosages"
[1] "Comparing accuracies:  0.716"