cpg_islands.cpp
上传用户:yhdzpy8989
上传日期:2007-06-13
资源大小:13604k
文件大小:10k
源码类别:

生物技术

开发平台:

C/C++

  1. /*
  2.  * ===========================================================================
  3.  * PRODUCTION $Log: cpg_islands.cpp,v $
  4.  * PRODUCTION Revision 1000.5  2004/06/01 20:54:52  gouriano
  5.  * PRODUCTION PRODUCTION: UPGRADED [GCC34_MSVC7] Dev-tree R1.20
  6.  * PRODUCTION
  7.  * ===========================================================================
  8.  */
  9. /*  $Id: cpg_islands.cpp,v 1000.5 2004/06/01 20:54:52 gouriano Exp $
  10. * ===========================================================================
  11. *
  12. *                            PUBLIC DOMAIN NOTICE
  13. *               National Center for Biotechnology Information
  14. *
  15. *  This software/database is a "United States Government Work" under the
  16. *  terms of the United States Copyright Act.  It was written as part of
  17. *  the author's official duties as a United States Government employee and
  18. *  thus cannot be copyrighted.  This software/database is freely available
  19. *  to the public for use. The National Library of Medicine and the U.S.
  20. *  Government have not placed any restriction on its use or reproduction.
  21. *
  22. *  Although all reasonable efforts have been taken to ensure the accuracy
  23. *  and reliability of the software and data, the NLM and the U.S.
  24. *  Government do not and cannot warrant the performance or results that
  25. *  may be obtained by using this software or data. The NLM and the U.S.
  26. *  Government disclaim all warranties, express or implied, including
  27. *  warranties of performance, merchantability or fitness for any particular
  28. *  purpose.
  29. *
  30. *  Please cite the author in any work or product based on this material.
  31. *
  32. * ===========================================================================
  33. *
  34. * Author: Philip Johnson
  35. *
  36. * File Description: cpg_island -- workbench algorithm plugin for
  37. * finding CpG islands
  38. *
  39. */
  40. #include <ncbi_pch.hpp>
  41. #include "cpg_islands.hpp"
  42. #include <algo/sequence/cpg.hpp>
  43. #include <gui/core/doc_manager.hpp>
  44. #include <gui/core/idocument.hpp>
  45. #include <gui/core/plugin_utils.hpp>
  46. #include <gui/core/version.hpp>
  47. #include <gui/plugin/PluginCommandSet.hpp>
  48. #include <gui/plugin/PluginInfo.hpp>
  49. #include <gui/plugin/PluginRequest.hpp>
  50. #include <gui/plugin/PluginValueConstraint.hpp>
  51. #include <gui/objutils/utils.hpp>
  52. #include <objects/seq/Seq_annot.hpp>
  53. #include <objects/seqfeat/Seq_feat.hpp>
  54. #include <objects/seqloc/Seq_interval.hpp>
  55. #include <serial/iterator.hpp>
  56. #include <objmgr/util/sequence.hpp>
  57. #include <objmgr/seq_vector.hpp>
  58. BEGIN_NCBI_SCOPE
  59. USING_SCOPE(objects);
  60. CAlgoPlugin_CpGIslands::~CAlgoPlugin_CpGIslands()
  61. {
  62. }
  63. // standard info boilerplate
  64. void CAlgoPlugin_CpGIslands::GetInfo(CPluginInfo& info)
  65. {
  66.     info.Reset();
  67.     // version info macro
  68.     info.SetInfo(CPluginVersion::eMajor, CPluginVersion::eMinor, 0,
  69.                  string(__DATE__) + " " + string(__TIME__),
  70.                  "CAlgoPlugin_CpGIslands", "Composition/CpG islands",
  71.                  "Scans for regions of unusually high CpG density.",
  72.                  "");
  73.     // command info
  74.     CPluginCommandSet& cmds = info.SetCommands();
  75.     CPluginCommand&    args = cmds.AddAlgoCommand(eAlgoCommand_run);
  76.     args.AddDefaultArgument("windowsize", "Size of sliding window",
  77.                             CPluginArg::eInteger, "200");
  78.     args.AddDefaultArgument("minlen", "Minimum length of island",
  79.                             CPluginArg::eInteger, "500");
  80.     args.AddDefaultArgument("gc", "Minimum %G+%C",
  81.                             CPluginArg::eDouble, "0.5");
  82.     args.AddDefaultArgument("cpg", "Minimum observed/expected CpG ratio",
  83.                             CPluginArg::eDouble, "0.6");
  84.     args.AddOptionalArgument("merge", "Merge adjacent islands within XXX",
  85.                              CPluginArg::eInteger);
  86.     args.AddArgument("locs", "Locations to scan", CSeq_loc::GetTypeInfo(),
  87.                      CPluginArg::TData::e_Array);
  88.     args.SetConstraint("locs",
  89.                        (*CPluginValueConstraint::CreateSeqMol(),
  90.                         CSeq_inst::eMol_na,
  91.                         CSeq_inst::eMol_dna,
  92.                         CSeq_inst::eMol_rna));
  93. }
  94. //-----------------------------------------------------------------------------
  95. // PRE : 
  96. // POST: 
  97. void CAlgoPlugin_CpGIslands::RunCommand(CPluginMessage& msg)
  98. {
  99.     const CPluginCommand& args = msg.GetRequest().GetCommand();
  100.     CPluginReply& reply = msg.SetReply();
  101.     if ( !m_Dialog.get() ) {
  102.         m_Dialog.reset(new CMultiColDlg());
  103.         m_Dialog->SetTitle("CpG Islands");
  104.         m_Dialog->SetColumn(0, "Sequence", FL_ALIGN_LEFT, 0.25f);
  105.         m_Dialog->SetColumn(1, "Location", FL_ALIGN_LEFT, 0.50f);
  106.         m_Dialog->SetColumn(2, "%GC",      FL_ALIGN_CENTER, 0.25f);
  107.     }
  108.     size_t sum_isles = 0;
  109.     size_t row = 0;
  110.     plugin_args::TLocList locs;
  111.     GetArgValue(args["locs"], locs);
  112.     ITERATE (plugin_args::TLocList, iter, locs) {
  113.         const CSeq_loc&  loc = *iter->second;
  114.         const IDocument& doc = *iter->first;
  115.         //first get bases
  116.         CBioseq_Handle handle = doc.GetScope().GetBioseqHandle(loc);
  117.         CSeqVector sv = handle
  118.             .GetSequenceView(loc, CBioseq_Handle::eViewConstructed,
  119.                              CBioseq_Handle::eCoding_Iupac);
  120.         string data;
  121.         sv.GetSeqData(0, sv.size(), data);
  122.         //find islands
  123.         CCpGIslands isles(data.data(), data.size(),
  124.                           args["windowsize"].AsInteger(),
  125.                           args["minlen"].AsInteger(),
  126.                           args["gc"].AsDouble(),
  127.                           args["cpg"].AsDouble());
  128.         //report islands
  129.         if (!isles.GetIsles().empty()) {                
  130.             sum_isles += isles.GetIsles().size();
  131.             string& id_str  = m_Dialog->SetCell(row, 0);
  132.             string& loc_str = m_Dialog->SetCell(row, 1);
  133.             const CSeq_id& best_id =
  134.                 sequence::GetId(handle, sequence::eGetId_Best);
  135.             id_str.erase();
  136.             best_id.GetLabel(&id_str);
  137.             loc_str = CPluginUtils::GetLabel(loc, &doc.GetScope());
  138.             CRef<CSeq_annot> sa(new CSeq_annot);
  139.             sa->AddName("CpG islands");
  140.             sa->AddTitle(string("CpG islands on ") + loc_str);
  141.             CNcbiOstrstream oss;
  142.             oss << "windowsize: " << args["windowsize"].AsInteger()
  143.                 << ", minlen: " << args["minlen"].AsInteger()
  144.                 << ", min %gc: " << args["gc"].AsDouble()
  145.                 << ", min observed/expected CpG: "
  146.                 << args["cpg"].AsDouble();
  147.             sa->AddComment(CNcbiOstrstreamToString(oss));
  148.             CSeq_annot::C_Data::TFtable &feats = sa->SetData().SetFtable();
  149.             ITERATE (CCpGIslands::TIsles, i, isles.GetIsles()) {
  150.                 CRef<CSeq_feat> feat(new CSeq_feat);
  151.                 CSeq_interval &seqInt = feat->SetLocation().SetInt();
  152.                 seqInt.SetFrom(i->m_Start);
  153.                 seqInt.SetTo(i->m_Stop);
  154.                 seqInt.SetId().Assign(sequence::GetId(loc));
  155.                 CRef<CSeq_loc> new_loc =
  156.                     CSeqUtils::RemapChildToParent(loc, feat->GetLocation());
  157.                 feat->SetLocation(*new_loc);
  158.                 size_t size = (i->m_Stop - i->m_Start);
  159.                 size_t pct_gc = (i->m_C + i->m_G) * 100 / size;
  160.                 string& gc_str  = m_Dialog->SetCell(row, 2);
  161.                 gc_str = NStr::IntToString(pct_gc);
  162.                 string& label = feat->SetData().SetRegion();
  163.                 label = "CpG island: " + gc_str;
  164.                 label += "% GC, ";
  165.                 label += NStr::IntToString(size);
  166.                 label += " bases";
  167.                 feats.push_back(feat);
  168.                 ++row;
  169.             }
  170.             reply.AddObject(doc, *sa);
  171.         }
  172.     }
  173.     if (reply.GetRaw().size() != 0) {
  174.         reply.AddAction(CPluginReplyAction::e_Add_to_document);
  175.     }
  176.     string label("Found ");
  177.     label += NStr::IntToString(sum_isles);
  178.     label += " isles on ";
  179.     label += NStr::IntToString(locs.size());
  180.     label += " locations:";
  181.     m_Dialog->SetLabel(label);
  182.     m_Dialog->SetRows(row);
  183.     m_Dialog->Show();
  184.     reply.SetStatus(eMessageStatus_success);
  185. }
  186. END_NCBI_SCOPE
  187. /*
  188. * ---------------------------------------------------------------------------
  189. * $Log: cpg_islands.cpp,v $
  190. * Revision 1000.5  2004/06/01 20:54:52  gouriano
  191. * PRODUCTION: UPGRADED [GCC34_MSVC7] Dev-tree R1.20
  192. *
  193. * Revision 1.20  2004/05/21 22:27:46  gorelenk
  194. * Added PCH ncbi_pch.hpp
  195. *
  196. * Revision 1.19  2004/05/03 13:05:42  dicuccio
  197. * gui/utils --> gui/objutils where needed
  198. *
  199. * Revision 1.18  2004/03/05 17:35:37  dicuccio
  200. * Use sequence::GetId() instead of CSeq_id::GetStringDescr()
  201. *
  202. * Revision 1.17  2004/01/27 18:37:36  dicuccio
  203. * Code clean-up.  Use standard names for plugins.  Removed unnecessary #includes
  204. *
  205. * Revision 1.16  2004/01/07 15:50:36  dicuccio
  206. * Adjusted for API change in CPluginUtils::GetLabel().  Standardized exception
  207. * reporting in algorithms.
  208. *
  209. * Revision 1.15  2003/12/12 20:07:08  johnson
  210. * accommodate MSVC 7 refactoring
  211. *
  212. * Revision 1.14  2003/11/24 15:45:25  dicuccio
  213. * Renamed CVersion to CPluginVersion
  214. *
  215. * Revision 1.13  2003/11/18 17:48:36  dicuccio
  216. * Added standard processing of return values
  217. *
  218. * Revision 1.12  2003/11/06 20:12:12  dicuccio
  219. * Cleaned up handling of USING_SCOPE - removed from all headers
  220. *
  221. * Revision 1.11  2003/11/04 17:49:22  dicuccio
  222. * Changed calling parameters for plugins - pass CPluginMessage instead of paired
  223. * CPluginCommand/CPluginReply
  224. *
  225. * Revision 1.10  2003/10/14 12:52:26  dicuccio
  226. * Fixed remapping of relative locations
  227. *
  228. * Revision 1.9  2003/10/07 13:47:00  dicuccio
  229. * Renamed CPluginURL* to CPluginValue*
  230. *
  231. * Revision 1.8  2003/09/04 14:05:24  dicuccio
  232. * Use IDocument instead of CDocument
  233. *
  234. * Revision 1.7  2003/09/03 14:46:53  rsmith
  235. * change namespace name from args to plugin_args to avoid clashes with variable names.
  236. *
  237. * Revision 1.6  2003/08/21 12:02:31  dicuccio
  238. * Added dialog box to display results.  Changed formatting of feature label
  239. *
  240. * Revision 1.5  2003/07/22 15:32:16  dicuccio
  241. * Changed to make use of new API in plugin_utils.hpp - GetArgValue()
  242. *
  243. * Revision 1.4  2003/07/21 19:32:53  dicuccio
  244. * Added constraints based on molecule type
  245. *
  246. * Revision 1.3  2003/07/19 13:53:35  ucko
  247. * Use CNcbiOstrstream rather than ostringstream, which GCC 2.9x lacks.
  248. *
  249. * Revision 1.2  2003/07/17 16:22:08  johnson
  250. * bug fix for msvc
  251. *
  252. * Revision 1.1  2003/07/16 17:19:40  johnson
  253. * Initial revision
  254. *
  255. * ===========================================================================
  256. */