Python db.Session類代碼示例

def xml(id):
    from hpf.hddb.db import Session, Family
    session = Session()
    family = session.query(Family).get(id)
    filename = "%i.xml" % family.id

    if runtime().opt(GZIP):
        import gzip
        filename = "%s.gz" % filename
        handle = gzip.open(filename,"w")
    else:
        handle = open(filename,"w")

    try:
        doc = FamilyFeatureBuilder(
            lambda: DefaultXMLGenerator(handle,pretty=True),
            lambda handler: StructureFeatureProvider(handler),
            lambda handler: ColumnFeatureProvider(handler),
            lambda handler: IeaFeatureProvider(handler),
            lambda handler: SelectionFeatureProvider(handler)
            )
        doc.buildDocument(family)
    finally:
        handle.close()
    session.close()

def populate_seqlength(sdict):
    from hpf.hddb.db import Session, Sequence
    session = Session()
    for seq_id in sdict.keys():
        seq_obj = session.query(Sequence).get(seq_id)
        sdict[seq_id].sequence_length = len(seq_obj.sequence)
    print "Getting sequence length for all sequences complete"

def astral_to_domain(experiment_id, threshold=0.5, dbstore=False):
    """
    Fetches all known-type domains for givein experiment ID. Computes astral overlap to
    domains, prints and optionally stores in hpf DB (table astral_domain_overlap)
    Will only store overlaps >= threshold parameter
    """
    from hpf.hddb.db import Session, push_to_db,  AstralDomainOverlap, Protein, Domain
    from hpf.structure_comparison.overlap import overlap
    from hpf.structure_comparison.astral_util import get_astrals, get_astral_startstop, parse_astral_chain

    # Create session and get all domains
    session = Session()
    domains = session.query(Domain).join(Protein).filter(Protein.experiment_key==experiment_id).filter(Domain.domain_type.in_(['psiblast','fold_recognition'])).all()
    print "Considering {0} domains to compute astral overlap for".format(len(domains))

    # For each domain, get representative astrals, calculate overlap, and store (optional)
    missing_astral = 0
    for domain in domains:
        domain_pdb_start = domain.region.parent_start
        domain_pdb_stop  = domain.region.parent_stop

        astrals = get_astrals(domain, session)
        if not astrals:
            #print "No astrals found for domain {0}".format(domain.id)
            missing_astral += 1
            continue

        for astral in astrals:
            try:
                (astral_start, astral_stop) = get_astral_startstop(astral)
                overlap_ratio = overlap(astral_start, astral_stop, domain_pdb_start, domain_pdb_stop)
            except ValueError:
                print "Negative overlap for domain {0} ({1}-{2}), Astral {3} (PDB {4}{5})".format(
                        domain.id, domain_pdb_start, domain_pdb_stop, astral.sid, astral.pdbid, astral.chain)
                print "Ignoring, moving to next astral.."
                continue
            except:
                print "Error calculating overlap for  domain {0} ({1}-{2}), Astral {3} (PBD {4}{5})".format(
                        domain.id, domain_pdb_start, domain_pdb_stop, astral.sid, astral.pdbid, astral.chain)
                raise

            if dbstore and overlap_ratio >= float(threshold):
                chain = parse_astral_chain(astral.chain)
                atod_dbo = AstralDomainOverlap(astral_id=astral.id, 
                                               astral_sid=astral.sid, 
                                               domain_id=domain.id, 
                                               astral_start=astral_start, 
                                               astral_stop=astral_stop, 
                                               domain_start=domain_pdb_start, 
                                               domain_stop=domain_pdb_stop, 
                                               pdb_id=astral.pdbid, 
                                               chain=chain, 
                                               overlap=overlap_ratio,
                                              )
                push_to_db(session, atod_dbo, exception_str="Error in pushing AstralDomainOverlap {0} to DB".format(atod_dbo), raise_on_duplicate=False)
            
            #print "Domain {0} ({1}-{2}), Astral {3} (PDB {4}{5}), Astral overlap {6}".format(domain.id, domain_pdb_start, domain_pdb_stop, astral.sid, astral.pdbid, astral.chain, overlap_ratio)

    print "Calculating astral to domain overlap for experiment {0} complete".format(experiment_id)
    print "{0} of {1} known-structure domains had no astral entries".format(missing_astral, len(domains))

class AlignmentToPDBMapper():
# Creates a mapper between a family alignment and a pdb structure
# NOTE: This is an abstraction of all of the confused and undocumented mappers below.
# Use this, because it is clear and because it works. KEEP IT SIMPLE.
# UNIT TESTS in pdb/tests/alignment_pdbmapper_test.py
# Access three public maps:
#   alignment_pdbseq_map  - alignment col => pdb seqres
#   pdbseq_pdbatom_map    - pdb seqres    => pdb atom
#   alignment_pdbatom_map - alignment col => pdb atom

    def __init__(self, family, protein, domain, debug=False):
        self.DEBUG = debug
        self.session = Session()
        self.family = family
        self.protein = protein
        self.domain  = domain
        self.alignment = self.family.alignments[0]
        if self.domain.sccs:
            self.pdbseqres = self.session.query(PDBSeqRes).filter_by(sequence_key=self.domain.parent_id[3:], chain=self.domain.sccs.chain).first()
        else:
            self.pdbseqres = self.session.query(PDBSeqRes).filter_by(sequence_key=self.domain.parent_id[3:]).first()
        self.pdbid = self.pdbseqres.pdb.pdbId+self.pdbseqres.chain

        # Create map: alignment column -> pdb seq res
        self.alignment_pdbseq_map = PDBDomainAlignmentMapper(self.protein, self.domain, self.alignment.alignment, self.pdbseqres, inverse=True)

        # Create map: pdb seq res -> pdb atom res
        self.pdbseq_pdbatom_map = PDBAtomSeqResMapper(self.pdbseqres, inverse=False)

        # Create map: alignment column -> pdb atom res
        self.alignment_pdbatom_map = AlignmentToPDBAtomMapper(self.alignment_pdbseq_map, self.pdbseq_pdbatom_map)

 def __init__(self, 
          protein, 
          domain, 
          alignment, 
          pdbseqres=None, 
          **kwargs):
     """
     @type protein: hpf.hddb.db.Protein
     @type domain: hpf.hddb.db.Domain
     @type alignment: Bio.Align.Generic
     @param pdbseqres: (Optional) Specify a PDB chain's sequence.
     """
     super(PDBDomainAlignmentMapper,self).__init__( protein, domain, alignment, **kwargs)
     self._chain = pdbseqres
     
     # If no pdbseqres given, fetch PDBSeqRes ORM object from db by sequence_key (a parsed domain.parent_id). 
     if self._chain==None:
         from hpf.hddb.db import Session,PDBSeqRes
         Session = Session()
         parent_id = int(self._domain.parent_id[3:])
         if domain.sccs:
             self._chain = Session.query(PDBSeqRes).filter_by(sequence_key=parent_id, chain=domain.sccs.chain).first()
         else:
             self._chain = Session.query(PDBSeqRes).filter(PDBSeqRes.sequence_key==parent_id).first()
     pdbid = self._chain.pdb.pdbId+self._chain.chain
     self._pdbid = pdbid;
     self._seed = self._seed_alignment()

    def setUp(self,):
        session = Session()
        # family_id  = 19187
        # protein_id = 1151960
        # domain_id  = 1302435
        family_id = 18883
        protein_id = 1063014
        domain_id = 1212855

        self.family = session.query(Family).get(family_id)
        self.protein = session.query(Protein).get(protein_id)
        self.domain = session.query(Domain).get(domain_id)

        self.ps_sites = get_possel_sites(self.family)
        (self.pdb_id, self.pdb_chain, self.struct) = get_domain_struct(self.domain)

        print "Testing on Family {0}, Protein {1}, Domain {2}, Structure {3}, PDB {4}{5}".format(
            family_id, protein_id, domain_id, self.struct.id, self.pdb_id, self.pdb_chain
        )
        print "Sites of +Sel for family: ", self.ps_sites

        self.pdb_struct = BioPDBStruct(self.pdb_id, self.pdb_chain, debug=True)
        self.cluster_id_str = "Family {0}, Protein {1}, Domain {2}, Structure {3}".format(
            family_id, protein_id, domain_id, self.struct.id
        )

 def __getitem__(self, key):
     from hpf.hddb.db import Session, SequenceAc
     amnh = self._oid_amnh[key]
     session = Session()
     try:
         ac = session.query(SequenceAc).filter(SequenceAc.ac==amnh).first()
         return ac.protein_key if ac else None
     finally:
         session.close()

    def setUp(self, ):
        session = Session()
        self.fam  = session.query(Family).get(19187)
        self.prot = session.query(Protein).get(1171456)
        self.dom  = session.query(Domain).get(1307995)
        self.alignment = self.fam.alignment.alignment
        self.psr  = session.query(PDBSeqRes).filter_by(sequence_key=self.dom.parent_id[3:]).first()

        self.ptpdba = ProteinToPDBAtom(self.prot,self.dom,self.alignment, pdbseqres=self.psr, inverse=False)
        self.pdbmap = AlignmentToPDBMapper(self.fam.id, self.prot.id, self.dom.id, debug=True)

def main():
    from hpf.hddb.db import Session, Sequence, Astral
    session = Session()
    astral_records = session.query(Astral)
    with open(args.outfile, 'w') as handle:
        for astral in astral_records:
            print "{0}\r".format(astral.sccs),
            handle.write(">{0}|{1}{2}|{3}\n".format(astral.sccs, astral.pdbid, astral.chain, astral.sequence_key))
            handle.write("{0}\n".format(astral.sequence.sequence))
    print "Exporting SCOP sequences to file {0} complete".format(args.outfile)

def init_dict():
# Creates a returns a dict of the form sequence id => Signaprint obj. Initially, signaprint object is empty
    from hpf.hddb.db import Session, Protein
    from sqlalchemy import distinct
    session = Session()
    mouse_seqs = session.query(distinct(Protein.sequence_key)).filter_by(experiment_key=1171)
    signaprint_dict = dict()
    for (seq,) in mouse_seqs:
        signaprint_dict[seq] = Signaprint(sequence_id=seq)
    print "Dict of {0} sequence => Signaprints created".format(len(signaprint_dict))
    return signaprint_dict

def ids_from_db():
# Retrieves a list of local (hpf db) sequence IDs for the sequences you want to calculate enrichment on
# These must be the same ID as in the ID field of pfam/interpro results

    # Get sequence IDs for all humanrna proteins (773 list)
    session = Session()
    proteins = session.query(Protein).filter_by(experiment_key=1177).all()
    protein_seqids = []
    for protein in proteins:
        protein_seqids.append(protein.sequence_key)
    return protein_seqids

    def setUp(self, ):
        self.family = 19187;
        self.protein = 1171456;
        self.domain  = 1307995;
        
        session = Session()
        self.family = session.query(Family).get(self.family)
        self.protein = session.query(Protein).get(self.protein)
        self.domain  = session.query(Domain).get(self.domain)

        self.pdbmap = AlignmentToPDBMapper(self.family, self.protein, self.domain, debug=True)

def index(records):
    """
    In-place rename's record id's from OID name to HPF sequence_key
    """
    from hpf.hddb.db import Session, SequenceAc
    mapping = oid_amnh()
    session = Session()
    for record in records:
        amnh = mapping[record.id]
        id = str(session.query(SequenceAc).filter(SequenceAc.ac==amnh).one().sequence_key)
        record.id = id

 def _make_seqfile(self, ):
     from hpf.hddb.db import Session, Sequence
     session=Session()
     sequence = session.query(Sequence).get(self.sequence_id)
     if not sequence:
         raise Exception("Getting sequence object from database failed")
     outhandle = open(self.sequence_file, 'w')
     outhandle.write(">hpf_seqid|{0}\n".format(sequence.id))
     outhandle.write("{0}\n".format(sequence.sequence))
     outhandle.close()
     session.close()

def tasks(experiment_id):
    print "Interpro driver::tasks:: Getting tasks (sequences) for experiment {0}".format(experiment_id)
    session = Session()
    proteins = session.query(Protein.sequence_key).filter_by(experiment_key=experiment_id)
    experiment_seqs = list()
    for protein in proteins:
        experiment_seqs.append(protein[0])
    if proteins.count() != len(experiment_seqs):
        raise Exception("Number of sequences extracted from experiment {0} does not match number of proteins".format(experiment_id))
    unique_seqs = list(set(experiment_seqs))
    print "{0} tasks (unique sequences) from {1} sequences retrieved".format(len(unique_seqs), len(experiment_seqs))
    return unique_seqs