Commit 5fcbe8b8 authored by Claus Jonathan Fritzemeier's avatar Claus Jonathan Fritzemeier
Browse files

groups for subsys integration

parent 82194f3c
......@@ -258,9 +258,9 @@ parseNotesReact <- function(notes) {
subSyst <- sub("SUBSYSTEM: *", "", fields_str[j])
subSyst <- sub("^S_", "", subSyst, perl = TRUE)
subSyst <- gsub("[_]+", " ", subSyst)
if (nchar(subSyst) == 0) {
subSyst <- "Exchange"
}
#if (nchar(subSyst) == 0) {
# subSyst <- "Exchange"
#}
#print(subSyst)
}
}
......@@ -646,6 +646,9 @@ for (i in 1 : numreact) {
}
# get subsystem properties from the sbml groups plugin
subSysGroups <- getSBMLGroupsList(sbmlmod)
# ---------------------------------------------------------------------------- #
# search for unused metabolites and unused reactions
......@@ -964,7 +967,9 @@ else {
sybil::gprRules(mod) <- rules
sybil::rxnGeneMat(mod) <- rxnGeneMat
#sybil::subSys(mod) <- subSys
sybil::subSys(mod) <- sybil:::.prepareSubSysMatrix(subSys, numreact)
if(is.null(subSysGroups)){
sybil::subSys(mod) <- sybil:::.prepareSubSysMatrix(subSys, numreact)
}
#sbml@gprRules <- rules
#sbml@genes <- genes
......@@ -1040,6 +1045,24 @@ if(newSybil)
if( !is.null(reactnotes) && length(reactnotes)==numreact )sybil::react_attr(mod)[['notes']]<-reactnotes
}
#------------------------------------------------------------------------------#
# subSys from groups #
#------------------------------------------------------------------------------#
if(!is.null(subSysGroups)){
subSysMat <- Matrix::Matrix(FALSE, nrow = numreact, ncol = length(subSysGroups), sparse = TRUE)
colnames(subSysMat) <- names(subSysGroups)
subSysMat <- Matrix::Matrix( do.call("cbind",
lapply(subSysGroups, function(x){
ids <- formatSBMLid(sub( "^R[_]+", "", x))
sybil::react_id(mod) %in% ids
})),
sparse = TRUE
)
sybil::subSys(mod) <- subSysMat
}
#------------------------------------------------------------------------------#
# Model Attr @Ardalan #
......
......@@ -222,6 +222,17 @@ getSBMLmodAnnotation <- function(sbmlm) {
#------------------------------------------------------------------------------#
getSBMLGroupsList <- function(sbmlm) {
modGroups <- .Call("getSBMLGroupsList", PACKAGE = "sybilSBML",
sbmlPointer(sbmlm)
)
return(modGroups)
}
#------------------------------------------------------------------------------#
getSBMLnumCompart <- function(sbmlm) {
num <- .Call("getSBMLnumCompart", PACKAGE = "sybilSBML",
......@@ -519,15 +530,21 @@ writeSBML<- function(morg=NULL,level=2,version=4,fbcLevel=0,filename="export.xml
} ####END newSybil attr
# Subsystem
if(is.null(newsubS) && !(modhasubSys) ) for ( i in 1:react_num(morg)) {newsubS[i]<- paste(names(which(subSys(morg)[i,])), collapse=", ")}
# format ids for sbml
newmet_id <- paste0("M_", (deformatSBMLid(met_id(morg))))
#newmet_id <- sub("\\[(\\w)\\]$", "_\\1", newmet_id) # append compartment id, if in postfix with square brkts
newreact_id <- paste0("R_", deformatSBMLid(react_id(morg)))
newmet_comp<-mod_compart(morg)[met_comp(morg)]
# Subsystem
if(is.null(newsubS) && !(modhasubSys) ) for ( i in 1:react_num(morg)) {newsubS[i]<- paste(names(which(subSys(morg)[i,])), collapse=", ")}
subSysGroups <- NULL
if(fbcLevel >= 2 && modhasubSys){
subSysGroups <- apply(subSys(morg), 2, function(x){
newreact_id[x]
})
names(subSysGroups) <- colnames(subSys(morg))
}
success <-.Call("exportSBML", PACKAGE = "sybilSBML",
as.integer(version),
......@@ -550,6 +567,7 @@ writeSBML<- function(morg=NULL,level=2,version=4,fbcLevel=0,filename="export.xml
as.numeric(uppbnd(morg)),
as.integer(obj_coef(morg)),
as.character(newsubS),
subSysGroups,
as.character(deformatGene(gpr(morg))),
as.numeric(shrinkMatrix(morg,j=1:react_num(morg))),
mod_notes,
......@@ -572,4 +590,4 @@ writeSBML<- function(morg=NULL,level=2,version=4,fbcLevel=0,filename="export.xml
}
else message(paste("Could not write file ",filename,"\n",sep=""), appendLF = FALSE);
return (success)
}
\ No newline at end of file
}
......@@ -42,6 +42,7 @@ static const R_CallMethodDef callMethods[] = {
{"getSBMLversion", (DL_FUNC) &getSBMLversion, 1},
{"validateDocument", (DL_FUNC) &validateDocument, 1},
{"getSBMLerrors", (DL_FUNC) &getSBMLerrors, 1},
{"getSBMLGroupsList", (DL_FUNC) &getSBMLGroupsList, 1},
{"getSBMLmodel", (DL_FUNC) &getSBMLmodel, 2},
{"getSBMLmodId", (DL_FUNC) &getSBMLmodId, 1},
{"getSBMLmodName", (DL_FUNC) &getSBMLmodName, 1},
......@@ -54,7 +55,7 @@ static const R_CallMethodDef callMethods[] = {
{"getSBMLCompartList", (DL_FUNC) &getSBMLCompartList, 1},
{"getSBMLSpeciesList", (DL_FUNC) &getSBMLSpeciesList, 1},
{"getSBMLReactionsList", (DL_FUNC) &getSBMLReactionsList, 1},
{"exportSBML", (DL_FUNC) &exportSBML, 33},
{"exportSBML", (DL_FUNC) &exportSBML, 34},
{"getSBMLFbcversion", (DL_FUNC) &getSBMLFbcversion, 1},
{NULL, NULL, 0}
};
......@@ -65,4 +66,4 @@ static const R_CallMethodDef callMethods[] = {
void R_init_sybilSBML(DllInfo *info) {
R_registerRoutines(info, NULL, callMethods, NULL, NULL);
R_useDynamicSymbols(info, FALSE);
}
\ No newline at end of file
}
......@@ -58,6 +58,17 @@ along with sybilSBML. If not, see <http://www.gnu.org/licenses/>.
#include <sbml/packages/fbc/sbml/FbcAnd.h>
#include <sbml/packages/fbc/sbml/FbcOr.h>
/*groups plugin*/
#include <sbml/packages/groups/common/GroupsExtensionTypes.h>
#include <sbml/packages/groups/extension/GroupsSBMLDocumentPlugin.h>
#include <sbml/packages/groups/extension/GroupsExtension.h>
#include <sbml/packages/groups/extension/GroupsModelPlugin.h>
#include <sbml/packages/groups/sbml/Group.h>
#include <sbml/packages/groups/sbml/ListOfGroups.h>
#include <sbml/packages/groups/sbml/ListOfMembers.h>
#include <sbml/packages/groups/sbml/Member.h>
static SEXP tagSBMLmodel;
static SEXP tagSBMLdocument;
......@@ -1161,6 +1172,41 @@ SEXP getSBMLSpeciesList(SEXP sbmlmod) {
}
SEXP getSBMLGroupsList(SEXP sbmlmod) {
GroupsModelPlugin_t * modelPlug = NULL;
modelPlug = (GroupsModelPlugin_t *) SBase_getPlugin((SBase_t *)(R_ExternalPtrAddr(sbmlmod)), "groups");
if(modelPlug != NULL){
int n = GroupsModelPlugin_getNumGroups(modelPlug);
SEXP rgroups = PROTECT(Rf_allocVector(VECSXP, n));
SEXP groupnames = PROTECT(Rf_allocVector(STRSXP, n));
for(int i=0; i<n; i++){
Group_t* group = GroupsModelPlugin_getGroup(modelPlug, i);
if(Group_isSetName(group) == 1){ // skip group if no name is set
SET_STRING_ELT(groupnames, i, Rf_mkChar(Group_getName(group)));
int m = Group_getNumMembers(group);
SEXP rmembers = PROTECT(Rf_allocVector(STRSXP, m));
for(int j=0; j < m; j++){
Member_t * member = Group_getMember(group, j);
char * memref = Member_getIdRef(member);
SET_STRING_ELT(rmembers, j, Rf_mkChar(memref));
}
SET_VECTOR_ELT(rgroups, i, rmembers);
UNPROTECT(1);
}else{
SET_VECTOR_ELT(rgroups, i, R_NilValue);
}
}
Rf_namesgets(rgroups, groupnames);
UNPROTECT(2);
return rgroups;
}else{
return R_NilValue;
}
return R_NilValue;
}
/* -------------------------------------------------------------------------- */
/* get list of reactions */
SEXP getSBMLReactionsList(SEXP sbmlmod) {
......@@ -1605,7 +1651,7 @@ void ParseModtoAnno (SBase_t* comp , char* Mannocopy)
SEXP exportSBML (SEXP version, SEXP level,SEXP FbcLevel, SEXP filename,SEXP sybil_max, SEXP mod_desc, SEXP mod_name, SEXP mod_compart, SEXP met_id, SEXP met_name, SEXP met_comp, SEXP met_form,SEXP met_charge, SEXP react_id, SEXP react_name, SEXP react_rev, SEXP lowbnd, SEXP uppbnd, SEXP obj_coef, SEXP subSys, SEXP gpr, SEXP SMatrix, SEXP mod_notes, SEXP mod_anno, SEXP com_notes , SEXP com_anno, SEXP met_notes, SEXP met_anno, SEXP met_bnd , SEXP react_notes, SEXP react_anno, SEXP ex_react,SEXP allgenes)
SEXP exportSBML (SEXP version, SEXP level, SEXP FbcLevel, SEXP filename, SEXP sybil_max, SEXP mod_desc, SEXP mod_name, SEXP mod_compart, SEXP met_id, SEXP met_name, SEXP met_comp, SEXP met_form, SEXP met_charge, SEXP react_id, SEXP react_name, SEXP react_rev, SEXP lowbnd, SEXP uppbnd, SEXP obj_coef, SEXP subSys, SEXP subSysGroups, SEXP gpr, SEXP SMatrix, SEXP mod_notes, SEXP mod_anno, SEXP com_notes , SEXP com_anno, SEXP met_notes, SEXP met_anno, SEXP met_bnd , SEXP react_notes, SEXP react_anno, SEXP ex_react, SEXP allgenes)
{
//Varaibles from R
const char* fname = CHAR(STRING_ELT(filename, 0));
......@@ -1615,6 +1661,7 @@ SEXP exportSBML (SEXP version, SEXP level,SEXP FbcLevel, SEXP filename,SEXP sybi
int SBMLlevel = INTEGER(level)[0];
int SBMLversion = INTEGER(version)[0];
int SBMLfbcversion = INTEGER(FbcLevel)[0];
int SBMLgroupsversion = 1;
double sybilmax = REAL(sybil_max)[0];
double sybilmin = sybilmax*(-1);
......@@ -1680,12 +1727,17 @@ SEXP exportSBML (SEXP version, SEXP level,SEXP FbcLevel, SEXP filename,SEXP sybi
SBMLExtension_t *sbmlext = SBMLExtensionRegistry_getExtension("fbc");
/* create the sbml namespaces object with fbc */
fbc = XMLNamespaces_create();
XMLNamespaces_t * fbc = XMLNamespaces_create();
XMLNamespaces_add(fbc, SBMLExtension_getURI(sbmlext, 3, 1, SBMLfbcversion), "fbc");
sbmlns = SBMLNamespaces_create(3, 1);
SBMLNamespaces_addNamespaces(sbmlns, fbc);
/* add groups extention */
SBMLExtension_t *sbmlgext = SBMLExtensionRegistry_getExtension("groups");
XMLNamespaces_t * groups = XMLNamespaces_create();
XMLNamespaces_add(groups, SBMLExtension_getURI(sbmlgext, SBMLlevel, SBMLversion, SBMLgroupsversion), "groups");
SBMLNamespaces_addNamespaces(sbmlns, groups);
/* create the document */
sbmlDoc = SBMLDocument_createWithSBMLNamespaces(sbmlns);
......@@ -2231,7 +2283,7 @@ SEXP exportSBML (SEXP version, SEXP level,SEXP FbcLevel, SEXP filename,SEXP sybi
const double *lower_bnd = REAL(lowbnd);
const double *upper_bnd = REAL(uppbnd);
char buf[20];
char buf[21]; // changed from 20 to 21 to avoid buffer overflow
// FBC1 FLUXBOUNDS
sprintf(buf, "LOWER_BOUND%d", i);
if (INTEGER(obj_coef)[i] != 1)
......@@ -2275,6 +2327,29 @@ SEXP exportSBML (SEXP version, SEXP level,SEXP FbcLevel, SEXP filename,SEXP sybi
}
}
/* add subsystem as groups if fbc is >= 2 */
if(SBMLfbcversion >= 2){
if(!Rf_isNull(subSysGroups)){
GroupsModelPlugin_t* groupsPlug = NULL;
groupsPlug = (GroupsModelPlugin_t*) SBase_getPlugin((SBase_t *)(model), "groups");
for(int i=0; i < Rf_length(subSysGroups); i++){
Group_t* newGroup = GroupsModelPlugin_createGroup(groupsPlug);
Group_setKindAsString(newGroup, "partonomy");
Group_setName(newGroup, CHAR(STRING_ELT(Rf_getAttrib(subSysGroups, R_NamesSymbol), i)));
SBase_setSBOTerm((SBase_t *) newGroup, 0000633);
for(int j=0; j < Rf_length(VECTOR_ELT(subSysGroups, i)); j++){
Member_t* newMember = Member_create(SBMLlevel, SBMLversion, SBMLgroupsversion);
Member_setIdRef(newMember, CHAR(STRING_ELT(VECTOR_ELT(subSysGroups, i), j)));
Group_addMember(newGroup, newMember);
}
//GroupsModelPlugin_addGroup(groupsPlug, newGroup);
}
}
}
// write SBML file
int result = writeSBML(sbmlDoc, fname);
SEXP out = R_NilValue;
......@@ -2287,4 +2362,4 @@ SEXP exportSBML (SEXP version, SEXP level,SEXP FbcLevel, SEXP filename,SEXP sybi
/* -------------------------------------------------------------------------- */
\ No newline at end of file
/* -------------------------------------------------------------------------- */
......@@ -105,8 +105,11 @@ SEXP getSBMLCompartList(SEXP sbmlmod);
/* get list of species (metabolites) */
SEXP getSBMLSpeciesList(SEXP sbmlmod);
/* get list of groups */
SEXP getSBMLGroupsList(SEXP sbmlmod);
/* get list of reactions */
SEXP getSBMLReactionsList(SEXP sbmlmod);
/* export Modelorg to SBML*/
SEXP exportSBML (SEXP version, SEXP level,SEXP FbcLevel, SEXP filename,SEXP sybil_max, SEXP mod_desc, SEXP mod_name, SEXP mod_compart, SEXP met_id, SEXP met_name, SEXP met_comp, SEXP met_form,SEXP met_charge, SEXP react_id, SEXP react_name, SEXP react_rev, SEXP lowbnd, SEXP uppbnd, SEXP obj_coef, SEXP subSys, SEXP gpr, SEXP SMatrix, SEXP mod_notes, SEXP mod_anno, SEXP com_notes , SEXP com_anno, SEXP met_notes, SEXP met_anno, SEXP met_bnd , SEXP react_notes, SEXP react_anno, SEXP ex_react, SEXP allgenes);
\ No newline at end of file
SEXP exportSBML (SEXP version, SEXP level, SEXP FbcLevel, SEXP filename, SEXP sybil_max, SEXP mod_desc, SEXP mod_name, SEXP mod_compart, SEXP met_id, SEXP met_name, SEXP met_comp, SEXP met_form, SEXP met_charge, SEXP react_id, SEXP react_name, SEXP react_rev, SEXP lowbnd, SEXP uppbnd, SEXP obj_coef, SEXP subSys, SEXP subSysGroups, SEXP gpr, SEXP SMatrix, SEXP mod_notes, SEXP mod_anno, SEXP com_notes , SEXP com_anno, SEXP met_notes, SEXP met_anno, SEXP met_bnd , SEXP react_notes, SEXP react_anno, SEXP ex_react, SEXP allgenes);
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment