Error in getModelFromHomology function

[SANDRA-ui-rgp]

Description of the issue:

An existing function, "getModelFromHomology," is not working. I am attempting to create a draft metabolic model of Plasmodium vivax. I have conducted BLAST searches on the sequences of P. berghei and P. vivax, which I have saved as "blast." "ipbeliver" is my reference model, with 'pberg' and 'pvp01' serving as the reference model and target model IDs, respectively.

Code:

modelpviv = getModelFromHomology({ipbeliver}, blast, 'pvp01',{'pberg'})

Error:
Standardizing grRules of template model with ID "pberg" ... done
The logical indices contain a true value outside of the array bounds.

Error in deleteUnusedGenes (line 45)
reducedModel.geneFrom=reducedModel.geneFrom(toKeep);

Error in getModelFromHomology (line 541)
draftModel=deleteUnusedGenes(draftModel,0);

System information

Installing from location C:\Users\Swapnam\AppData\Roaming\MathWorks\MATLAB Add-Ons\Collections\RAVEN Toolbox
Checking RAVEN release 2.8.6
Latest RAVEN release available 2.8.7
Instructions on how to upgrade here
Checking MATLAB release 2023b
Checking system architecture win64
Set RAVEN in MATLAB path Pass
Save MATLAB path Pass

[edkerk]

To troubleshoot this, please do the following:

• Run open getModelFromHomology
• Add %-symbols in the beginning of lines 541 to 544 (commenting out the commands on those lines)
• Save getModelFromHomology
• Rerun your code
• Once finished, run save('troubleshoot.mat','modelpviv')
• Post the generated mat file in a comment in this issue (you can remove it again later)

We can then try to replicate where the error is encountered.

[SANDRA-ui-rgp]

Thank you for your response, here is the mat file.
troubleshoot.zip

[edkerk]

I'm not sure why, but the geneFrom field (which is inherited from the template model) has less entries than the number of genes in the new model. Regardless, the geneFrom has lost its meaning after running getModelFromHomology, as the original genes are replaced by their orthologs. It is therefore most sensible to remove this field.

I will modify getModelFromHomology to remove the geneFrom field automatically, but as a workaround you should run this line of code before running getModelFromHomology:

ipbeliver = rmfield(ipbeliver,'geneFrom');