% ChemmineR
% Yiqun Cao, Tyler Backman, Kevin Horan, Thomas Girke
% Last update: `r format(Sys.time(),"%B %d, %Y")`
```{r echo=FALSE, results='hide',message=FALSE,warning=FALSE,error=FALSE}
library(ChemmineR)
library(knitcitations)
RefManageR:::BibOptions(style="markdown",no.print.fields=c("doi"))
sapply(read.bibtex(file="references.bib"),record_as_cited)
```
ChemmineR
=========
Yiqun Cao, Tyler Backman, Kevin Horan, Thomas Girke
Last update: `r format(Sys.time(),"%B %d, %Y")`
Email contact: thomas.girke@ucr.edu
Introduction
============
*ChemmineR* is a cheminformatics package for analyzing
drug-like small molecule data in R. Its latest version contains
functions for efficient processing of large numbers of small molecules,
physicochemical/structural property predictions, structural similarity
searching, classification and clustering of compound libraries with a
wide spectrum of algorithms.
In addition, *ChemmineR* offers visualization functions
for compound clustering results and chemical structures. The integration
of chemoinformatic tools with the R programming environment has many
advantages, such as easy access to a wide spectrum of statistical
methods, machine learning algorithms and graphic utilities. The first
version of this package was published in
`r citet(c("Cao_2008"))` . Since then many additional
utilities and add-on packages have been added to the environment and
many more are under development for future releases (Figure 2;
`r citep(c("Backman_2011","Wang_2013"))`).
- Improved SMILES support via new *SMIset* object class
and SMILES import/export functions
- Integration of a subset of OpenBabel functionalities via new
*ChemmineOB* add-on package `r citep(c("Cao_2008"))`
- Streaming functionality for processing millions of molecules on a
laptop
- Mismatch tolerant maximum common substructure (MCS) search algorithm
- Fast and memory efficient fingerprint search support using atom pair
or PubChem fingerprints
Getting Started
===============
Installation
------------
The R software for running ChemmineR can be downloaded from CRAN
(). The ChemmineR package can be
installed from R using the *bioLite* install command.
```{r eval=FALSE}
source("http://bioconductor.org/biocLite.R") # Sources the biocLite.R installation script.
biocLite("ChemmineR") # Installs the package.
```
Loading the Package and Documentation
-------------------------------------
```{r eval=TRUE, tidy=FALSE}
library("ChemmineR") # Loads the package
```
```{r eval=FALSE, tidy=FALSE}
library(help="ChemmineR") # Lists all functions and classes
vignette("ChemmineR") # Opens this PDF manual from R
```
Five Minute Tutorial
--------------------
The following code gives an overview of the most important
functionalities provided by *ChemmineR*. Copy and paste
of the commands into the R console will demonstrate their utilities.
Create Instances of *SDFset* class:
```{r eval=TRUE, tidy=FALSE}
data(sdfsample)
sdfset <- sdfsample
sdfset # Returns summary of SDFset
sdfset[1:4] # Subsetting of object
sdfset[[1]] # Returns summarized content of one SDF
```{r eval=FALSE, tidy=FALSE}
view(sdfset[1:4]) # Returns summarized content of many SDFs, not printed here
as(sdfset[1:4], "list") # Returns complete content of many SDFs, not printed here
```
An *SDFset* is created during the import of an SD file:
```{r eval=FALSE, tidy=FALSE}
sdfset <- read.SDFset("http://faculty.ucr.edu/ tgirke/Documents/R_BioCond/Samples/sdfsample.sdf")
```
Miscellaneous accessor methods for *SDFset* container:
```{r eval=FALSE, tidy=FALSE}
header(sdfset[1:4]) # Not printed here
```
```{r eval=TRUE, tidy=FALSE}
header(sdfset[[1]])
```
```{r eval=FALSE, tidy=FALSE}
atomblock(sdfset[1:4]) # Not printed here
```
```{r eval=TRUE, tidy=FALSE}
atomblock(sdfset[[1]])[1:4,]
```
```{r eval=FALSE, tidy=FALSE}
bondblock(sdfset[1:4]) # Not printed here
```
```{r eval=TRUE, tidy=FALSE}
bondblock(sdfset[[1]])[1:4,]
```
```{r eval=FALSE, tidy=FALSE}
datablock(sdfset[1:4]) # Not printed here
```
```{r eval=TRUE, tidy=FALSE}
datablock(sdfset[[1]])[1:4]
```
Assigning compound IDs and keeping them unique:
```{r eval=TRUE, tidy=FALSE}
cid(sdfset)[1:4] # Returns IDs from SDFset object
sdfid(sdfset)[1:4] # Returns IDs from SD file header block
unique_ids <- makeUnique(sdfid(sdfset))
cid(sdfset) <- unique_ids
```
Converting the data blocks in an *SDFset* to a matrix:
```{r eval=TRUE, tidy=FALSE}
blockmatrix <- datablock2ma(datablocklist=datablock(sdfset)) # Converts data block to matrix
numchar <- splitNumChar(blockmatrix=blockmatrix) # Splits to numeric and character matrix
numchar[[1]][1:2,1:2] # Slice of numeric matrix
numchar[[2]][1:2,10:11] # Slice of character matrix
```
Compute atom frequency matrix, molecular weight and formula:
```{r eval=TRUE, tidy=FALSE}
propma <- data.frame(MF=MF(sdfset), MW=MW(sdfset), atomcountMA(sdfset))
propma[1:4, ]
```
Assign matrix data to data block:
```{r eval=TRUE, tidy=FALSE}
datablock(sdfset) <- propma
datablock(sdfset[1])
```
String searching in *SDFset* ():
```{r eval=FALSE, tidy=FALSE}
grepSDFset("650001", sdfset, field="datablock", mode="subset") # Returns summary view of matches. Not printed here.
```
```{r eval=TRUE, tidy=FALSE}
grepSDFset("650001", sdfset, field="datablock", mode="index")
```
Export SDFset to SD file:
```{r eval=FALSE, tidy=FALSE}
write.SDF(sdfset[1:4], file="sub.sdf", sig=TRUE)
```
Plot molecule structure of one or many SDFs:
```{r plotstruct, eval=TRUE, tidy=FALSE}
plot(sdfset[1:4], print=FALSE) # Plots structures to R graphics device
```
```{r eval=FALSE, tidy=FALSE}
sdf.visualize(sdfset[1:4]) # Compound viewing in web browser
```
Structure similarity searching and clustering:
```{r eval=FALSE, tidy=FALSE}
apset <- sdf2ap(sdfset) # Generate atom pair descriptor database for searching
```
```{r eval=TRUE, tidy=FALSE}
data(apset) # Load sample apset data provided by library.
cmp.search(apset, apset[1], type=3, cutoff = 0.3, quiet=TRUE) # Search apset database with single compound.
cmp.cluster(db=apset, cutoff = c(0.65, 0.5), quiet=TRUE)[1:4,] # Binning clustering using variable similarity cutoffs.
```
OpenBabel Functions
===================
*ChemmineR* integrates now a subset of cheminformatics
functionalities implemented in the OpenBabel C++ library
`r citep(c("greycite13432", "Cao_2008" ))` . These
utilities can be accessed by installing the *ChemmineOB*
package and the OpenBabel software itself. *ChemmineR*
will automatically detect the availability of
*ChemmineOB* and make use of the additional utilities.
The following lists the functions and methods that make use of
OpenBabel. References are included to locate the sections in the manual
where the utility and usage of these functions is described.
*Structure format interconversions* (see Section [Format Inter-Conversions](#format-interconversions))
- *smiles2sdf*: converts from SMILES to SDF object
- *sdf2smiles*: converts from SDF to SMILES object
- *convertFormat*: converts strings between two formats
- *convertFormatFile*: converts files between two formats. This function can be used to enable ChemmineR to read in any
format supported by Open Babel. For example, if you had an SML file you could do:
```{r eval=FALSE, tidy=FALSE}
convertFormatFile("SML","SDF","mycompound.sml","mycompound.sdf")
sdfset=read.SDFset("mycompound.sdf")
```
*propOB*: generates several compound properties. See the man page for a current list of properties computed.
```{r eval=TRUE, tidy=FALSE}
propOB(sdfset[1])
```
*fingerprintOB*: generates fingerprints for compounds. The fingerprint name can be anything supported by OpenBabel. See the man page
for a current list.
```{r eval=TRUE, tidy=FALSE}
fingerprintOB(sdfset,"FP2")
```
*smartsSearchOB*: find matches of SMARTS patterns in compounds
```{r eval=TRUE, tidy=FALSE}
#count rotable bonds
smartsSearchOB(sdfset[1:5],"[!$(*#*)&!D1]-!@[!$(*#*)&!D1]",uniqueMatches=FALSE)
```
*exactMassOB*: Compute the monoisotopic (exact) mass of a set of compounds
```{r eval=TRUE, tidy=FALSE}
exactMassOB(sdfset[1:5])
```
*regenerateCoords*: Re-compute the 2D coordinates of a compound using Open Babel. This can sometimes
improve the quality of the compounds plot. See also the *regenCoords* option of the plot function.
```{r eval=TRUE, tidy=FALSE}
sdfset2 = regenerateCoords(sdfset[1:5])
plot(sdfset[1], regenCoords=TRUE,print=FALSE)
```
*generate3DCoords*: Generate 3D coordinates for compounds with only 2D coordinates.
```{r eval=FALSE, tidy=FALSE}
sdf3D = generate3DCoords(sdfset[1])
```
*canonicalize*: Compute a canonicalized atom numbering. This allows compounds with the same molecular
structure but different atom numberings to be compared properly.
```{r eval=TRUE, tidy=FALSE}
canonicalSdf= canonicalize(sdfset[1])
```
*canonicalNumbering*: Return a mapping from the original atom numbering to the
canonical atom number.
```{r eval=TRUE, tidy=FALSE}
mapping = canonicalNumbering(sdfset[1])
```
Overview of Classes and Functions
=================================
The following list gives an overview of the most important S4 classes,
methods and functions available in the ChemmineR package. The help
documents of the package provide much more detailed information on each
utility. The standard R help documents for these utilities can be
accessed with this syntax: *?function\_name* (*e.g.*
*?cid*) and *?class\_name-class* (*e.g*.
*?"SDFset-class"*).
Molecular Structure Data
------------------------
Classes
- *SDFstr*: intermediate string class to facilitate SD
file import; not important for end user
- *SDF*: container for single molecule imported from an
SD file
- *SDFset*: container for many SDF objects; most
important structure container for end user
- *SMI*: container for a single SMILES string
- *SMIset*: container for many SMILES strings
Functions/Methods (mainly for *SDFset* container,
*SMIset* should be coerced with
*smiles2sdf* to *SDFset*)
- Accessor methods for *SDF/SDFset*
- Object slots: *cid, header, atomblock, bondblock,
datablock (sdfid, datablocktag)*
- Summary of *SDFset*: *view*
- Matrix conversion of data block: *datablock2ma,
splitNumChar*
- String search in SDFset: *grepSDFset*
- Coerce one class to another
- Standard syntax *as(..., "...")* works in most
cases. For details see R help with
*?"SDFset-class"*.
- Utilities
- Atom frequencies: *atomcountMA, atomcount*
- Molecular weight: *MW*
- Molecular formula: *MF*
- ...
- Compound structure depictions
- R graphics device: *plot, plotStruc*
- Online: *cmp.visualize*
Structure Descriptor Data
-------------------------
Classes
- *AP*: container for atom pair descriptors of a single
molecule
- *APset*: container for many AP objects; most
important structure descriptor container for end user
- *FP*: container for fingerprint of a single molecule
- *FPset*: container for fingerprints of many
molecules, most important structure descriptor container for end
user
Functions/Methods
- Create *AP/APset* instances
- From *SDFset*: *sdf2ap*
- From SD file: *cmp.parse*
- Summary of *AP/APset*: *view,
db.explain*
- Accessor methods for AP/APset
- Object slots: *ap, cid*
- Coerce one class to another
- Standard syntax *as(..., "...")* works in most
cases. For details see R help with
*?"APset-class"*.
- Structure Similarity comparisons and Searching
- Compute pairwise similarities : *cmp.similarity,
fpSim*
- Search APset database: *cmp.search, fpSim*
- AP-based Structure Similarity Clustering
- Single-linkage binning clustering: *cmp.cluster*
- Visualize clustering result with MDS:
*cluster.visualize*
- Size distribution of clusters: *cluster.sizestat*
- Folding
- Fold a descriptor with *fold*
- Query the number of times a descriptor has been folded:
*foldCount*
- Query the number of bits in a descriptor: *numBits*
Import of Compounds
===================
SDF Import
----------
The following gives an overview of the most important import/export
functionalities for small molecules provided by
*ChemmineR*. The given example creates an instance of the
*SDFset* class using as sample data set the first 100
compounds from this PubChem SD file (SDF):
Compound\_00650001\_00675000.sdf.gz
().
SDFs can be imported with the *read.SDFset* function:
```{r eval=FALSE, tidy=FALSE}
sdfset <- read.SDFset("http://faculty.ucr.edu/ tgirke/Documents/R_BioCond/Samples/sdfsample.sdf")
```
```{r eval=TRUE, tidy=FALSE}
data(sdfsample) # Loads the same SDFset provided by the library
sdfset <- sdfsample
valid <- validSDF(sdfset) # Identifies invalid SDFs in SDFset objects
sdfset <- sdfset[valid] # Removes invalid SDFs, if there are any
```
Import SD file into *SDFstr* container:
```{r eval=FALSE, tidy=FALSE}
sdfstr <- read.SDFstr("http://faculty.ucr.edu/ tgirke/Documents/R_BioCond/Samples/sdfsample.sdf")
```
Create
*SDFset* from *SDFstr* class:
```{r eval=TRUE, tidy=FALSE}
sdfstr <- as(sdfset, "SDFstr")
sdfstr
as(sdfstr, "SDFset")
```
SMILES Import
-------------
The *read.SMIset* function imports one or many molecules
from a SMILES file and stores them in a *SMIset*
container. The input file is expected to contain one SMILES string per
row with tab-separated compound identifiers at the end of each line. The
compound identifiers are optional.
Create sample SMILES file and then import it:
```{r eval=FALSE, tidy=FALSE}
data(smisample); smiset <- smisample
write.SMI(smiset[1:4], file="sub.smi")
smiset <- read.SMIset("sub.smi")
```
Inspect content of *SMIset*:
```{r eval=TRUE, tidy=FALSE}
data(smisample) # Loads the same SMIset provided by the library
smiset <- smisample
smiset
view(smiset[1:2])
```
Accessor functions:
```{r eval=TRUE, tidy=FALSE}
cid(smiset[1:4])
smi <- as.character(smiset[1:2])
```
Create *SMIset* from named character vector:
```{r eval=TRUE, tidy=FALSE}
as(smi, "SMIset")
```
Export of Compounds
===================
SDF Export
----------
Write objects of classes *SDFset/SDFstr/SDF* to SD file:
```{r eval=FALSE, tidy=FALSE}
write.SDF(sdfset[1:4], file="sub.sdf")
```
Writing customized *SDFset* to file containing
*ChemmineR* signature, IDs from *SDFset*
and no data block:
```{r eval=FALSE, tidy=FALSE}
write.SDF(sdfset[1:4], file="sub.sdf", sig=TRUE, cid=TRUE, db=NULL)
```
Example for injecting a custom matrix/data frame into the data block of
an *SDFset* and then writing it to an SD file:
```{r eval=FALSE, tidy=FALSE}
props <- data.frame(MF=MF(sdfset), MW=MW(sdfset), atomcountMA(sdfset))
datablock(sdfset) <- props
write.SDF(sdfset[1:4], file="sub.sdf", sig=TRUE, cid=TRUE)
```
Indirect export via *SDFstr* object:
```{r eval=FALSE, tidy=FALSE}
sdf2str(sdf=sdfset[[1]], sig=TRUE, cid=TRUE) # Uses default components
sdf2str(sdf=sdfset[[1]], head=letters[1:4], db=NULL) # Uses custom components for header and data block
```
Write *SDF*, *SDFset* or
*SDFstr* classes to file:
```{r eval=FALSE, tidy=FALSE}
write.SDF(sdfset[1:4], file="sub.sdf", sig=TRUE, cid=TRUE, db=NULL)
write.SDF(sdfstr[1:4], file="sub.sdf")
cat(unlist(as(sdfstr[1:4], "list")), file="sub.sdf", sep="")
```
SMILES Export
-------------
Write objects of class *SMIset* to SMILES file with and
without compound identifiers:
```{r eval=FALSE, tidy=FALSE}
data(smisample); smiset <- smisample # Sample data set
write.SMI(smiset[1:4], file="sub.smi", cid=TRUE) write.SMI(smiset[1:4], file="sub.smi", cid=FALSE)
```
Format Interconversions
=======================
The *sdf2smiles* and *smiles2sdf*
functions provide format interconversion between SMILES strings
(Simplified Molecular Input Line Entry Specification) and
*SDFset* containers.
Convert an *SDFset* container to a SMILES
*character* string:
```{r sdf2smiles, eval=FALSE, tidy=FALSE}
data(sdfsample);
sdfset <- sdfsample[1]
smiles <- sdf2smiles(sdfset)
smiles
```
Convert a SMILES *character* string to an
*SDFset* container:
```{r smiles2sdf, eval=FALSE, tidy=FALSE}
sdf <- smiles2sdf("CC(=O)OC1=CC=CC=C1C(=O)O")
view(sdf)
```
When the *ChemineOB* package is installed these
conversions are performed with the OpenBabel Open Source Chemistry
Toolbox. Otherwise the functions will fall back to using the ChemMine
Tools web service for this operation. The latter will require internet
connectivity and is limited to only the first compound given.
*ChemmineOB* provides access to the compound format
conversion functions of OpenBabel. Currently, over 160 formats are
supported by OpenBabel. The functions *convertFormat* and
*convertFormatFile* can be used to convert files or
strings between any two formats supported by OpenBabel. For example, to
convert a SMILES string to an SDF string, one can use the
*convertFormat* function.
```{r eval=FALSE, tidy=FALSE}
sdfStr <- convertFormat("SMI","SDF","CC(=O)OC1=CC=CC=C1C(=O)O_name")
```
This will return the given compound as an SDF formatted string. 2D
coordinates are also computed and included in the resulting SDF string.
To convert a file with compounds encoded in one format to another
format, the *convertFormatFile* function can be used
instead.
```{r eval=FALSE, tidy=FALSE}
convertFormatFile("SMI","SDF","test.smiles","test.sdf")
```
To see the whole list of file formats supported by OpenBabel, one can
run from the command-line "obabel -L formats".
Splitting SD Files
==================
The following *write.SDFsplit* function allows to split
SD Files into any number of smaller SD Files. This can become important
when working with very big SD Files. Users should note that this
function can output many files, thus one should run it in a dedicated
directory!
Create sample SD File with 100 molecules:
```{r eval=FALSE, tidy=FALSE}
write.SDF(sdfset, "test.sdf")
```
Read in sample SD File. Note: reading file into SDFstr is much faster
than into SDFset:
```{r eval=FALSE, tidy=FALSE}
sdfstr <- read.SDFstr("test.sdf")
```
Run export on *SDFstr* object:
```{r eval=FALSE, tidy=FALSE}
write.SDFsplit(x=sdfstr, filetag="myfile", nmol=10) # 'nmol' defines the number of molecules to write to each file
```
Run export on *SDFset* object:
```{r eval=FALSE, tidy=FALSE}
write.SDFsplit(x=sdfset, filetag="myfile", nmol=10)
```
Streaming Through Large SD Files
================================
The *sdfStream* function allows to stream through SD
Files with millions of molecules without consuming much memory. During
this process any set of descriptors, supported by
*ChemmineR*, can be computed (*e.g.* atom pairs,
molecular properties, etc.), as long as they can be returned in tabular
format. In addition to descriptor values, the function returns a line
index that gives the start and end positions of each molecule in the
source SD File. This line index can be used by the downstream
*read.SDFindex* function to retrieve specific molecules
of interest from the source SD File without reading the entire file into
R. The following outlines the typical workflow of this streaming
functionality in *ChemmineR*.
Create sample SD File with 100 molecules:
```{r eval=FALSE, tidy=FALSE}
write.SDF(sdfset, "test.sdf")
```
Define descriptor set in a simple function:
```{r eval=FALSE, tidy=FALSE}
desc <- function(sdfset)
cbind(SDFID=sdfid(sdfset),
# datablock2ma(datablocklist=datablock(sdfset)),
MW=MW(sdfset),
groups(sdfset), APFP=desc2fp(x=sdf2ap(sdfset), descnames=1024,
type="character"), AP=sdf2ap(sdfset, type="character"), rings(sdfset,
type="count", upper=6, arom=TRUE) )
```
Run *sdfStream* with *desc* function and
write results to a file called *matrix.xls*:
```{r eval=FALSE, tidy=FALSE}
sdfStream(input="test.sdf", output="matrix.xls", fct=desc, Nlines=1000) # 'Nlines': number of lines to read from input SD File at a time
```
One can also start reading from a specific line number in the SD file.
The following example starts at line number 950. This is useful for
restarting and debugging the process. With *append=TRUE*
the result can be appended to an existing file.
```{r eval=FALSE, tidy=FALSE}
sdfStream(input="test.sdf", output="matrix2.xls", append=FALSE, fct=desc, Nlines=1000, startline=950)
```
Select molecules meeting certain property criteria from SD File using
line index generated by previous *sdfStream* step:
```{r eval=FALSE, tidy=FALSE}
indexDF <- read.delim("matrix.xls", row.names=1)[,1:4]
indexDFsub <- indexDF[indexDF$MW < 400, ] # Selects molecules with MW < 400
sdfset <- read.SDFindex(file="test.sdf", index=indexDFsub, type="SDFset") # Collects results in 'SDFset' container
```
Write results directly to SD file without storing larger numbers of
molecules in memory:
```{r eval=FALSE, tidy=FALSE}
read.SDFindex(file="test.sdf", index=indexDFsub, type="file",
outfile="sub.sdf")
```
Read AP/APFP strings from file into *APset* or
*FP* object:
```{r eval=FALSE, tidy=FALSE}
apset <- read.AP(x="matrix.xls", type="ap", colid="AP")
apfp <- read.AP(x="matrix.xls", type="fp", colid="APFP")
```
Alternatively, one can provide the AP/APFP strings in a named character
vector:
```{r eval=FALSE, tidy=FALSE}
apset <- read.AP(x=sdf2ap(sdfset[1:20], type="character"), type="ap")
fpchar <- desc2fp(sdf2ap(sdfset[1:20]), descnames=1024, type="character")
fpset <- as(fpchar, "FPset")
```
Storing Compounds in an SQL Database
====================================
As an alternative to sdfStream, there is now also an option to store
data in an SQL database, which then allows for fast queries and compound
retrieval. The default database
is SQLite, but any other SQL database should work with some minor
modifications to the table definitions, which are stored in
schema/compounds.SQLite under the ChemmineR package directory. Compounds
are stored in their entirety in the databases so there is no need to
keep any original data files.
Users can define their own set of compound features to compute and store
when loading new compounds. Each of these features will be stored in its
own, indexed table. Searches can then be performed using these features
to quickly find specific compounds. Compounds can always be retrieved
quickly because of the database index, no need to scan a large compound
file. In addition to user defined features, descriptors can also be
computed and stored for each compound.
A new database can be created with the *initDb* function.
This takes either an existing database connection, or a filename. If a
filename is given then an SQLite database connection is created. It then
ensures that the required tables exist and creates them if not. The
connection object is then returned. This function can be called safely
on the same connection or database many times and will not delete any
data.
Loading Data
------------
The functions *loadSdf* and *loadSmiles* can be used to load
compound data from either a file (both) or an *SDFset* (*loadSdf*
only). The *fct* parameter should be a function to
extract features from the data. It will be handed an
*SDFset* generated from the data being loaded. This may
be done in batches, so there is no guarantee that the given SDFSset will
contain the whole dataset. This function should return a data frame with
a column for each feature and a row for each compound given. The order
of the final data frame should be the same as that of the
*SDFset*. The column names will become the feature names.
Each of these features will become a new, indexed, table in the database
which can be used later to search for compounds.
The *descriptors* parameter can be a function which
computes descriptors. This function will also be given an
*SDFset* object, which may be done in batches. It should
return a data frame with the following two columns: "descriptor" and
"descriptor\_type". The "descriptor" column should contain a string
representation of the descriptor, and "descriptor\_type" is the type of
the descriptor. Our convention for atom pair is "ap" and "fp" for finger
print. The order should also be maintained.
When the data has been loaded, *loadSdf* will return the
compound id numbers of each compound loaded. These compound id numbers
are computed by the database and are not extracted from the compound
data itself. They can be used to quickly retrieve compounds later.
New features can also be added using this function. However, all
compounds must have all features so if new features are added to a new
set of compounds, all existing features must be computable by the
*fct* function given. If new features are detected, all
existing compounds will be run through *fct* in order to
compute the new features for them as well.
For example, if dataset X is loaded with features F1 and F2, and then at
a later time we load dataset Y with new feature F3, the
*fct* function used to load dataset Y must compute and
return features F1, F2, and F3. *loadSdf* will call
*fct* with both datasets X and Y so that all features are
available for all compounds. If any features are missing an error will
be raised. If just new features are being added, but no new compounds,
use the *addNewFeatures* function.
In this example, we create a new database called "test.db" and load it
with data from an *SDFset*. We also define
*fct* to compute the molecular weight, "MW", and the
number of rings and aromatic rings. The rings function actually returns
a data frame with columns "RINGS" and "AROMATIC", which will be merged
into the data frame being created which will also contain the "MW"
column. These will be the names used for these features and must be used
when searching with them. Finally, the new compound ids are returned and
stored in the "ids" variable.
```{r eval=TRUE, tidy=FALSE}
data(sdfsample)
#create and initialize a new SQLite database
conn <- initDb("test.db")
# load data and compute 3 features: molecular weight, with the MW function,
# and counts for RINGS and AROMATIC, as computed by rings, which
# returns a data frame itself.
ids<-loadSdf(conn,sdfsample, function(sdfset)
data.frame(rings(sdfset,type="count",upper=6, arom=TRUE),propOB(sdfset)) )
#list features in the database:
print(listFeatures(conn))
```
Updates
-------
By default the *loadSdf* / *loadSmiles* functions will detect duplicate
compound entries and only insert one of them. This means it is safe
to run these functions on the same data set several times and you
won't end up with duplicates. This allows the functions to be re-run
in the event that a previous run on a dataset does not complete.
Duplicate compounds are detected by compouting the MD5 checksum on
the textual representation of it.
It can also update existing compounds with new versions of the same
compound. To enable this, set *updateByName* to true. It will then
consider two compounds with the same name to be the same, even if the
definition is different. Then, if the name of a compound exists in
the database and it is trying to insert another compound with the
same name, it will overwrite the existing compound. It will also drop
and re-compute all associated descriptors and features for the new
compound (assuming the required functions for descriptor and feature
computation are available at the time the update is performed).
Duplicate Descriptors
---------------------
It is often the case when loading a large set of compounds that
several compounds will produce the same descriptor. *ChemmineR*
detects this case and only stores one copy of the descriptor for
every compound it is for. This feature saves some space and some
time for processes that need to be applied to every descriptor.
It also highlights a new problem. If you have a descriptor in hand
and you want to find a single compound to represent it, which
compound should be used if the descriptor was produced from multiple
compounds? To address this problem, *ChemmineR* allows you to set
priority values for each compound-descriptor mapping. Then, in
contexts where a single compound is required, the highest priority
compound will be chosen. Highest priority corresponds to the lowest
numerical value. So mapping with priority 0 would be used first.
To set these priorities there is the function *setPriorities*.
It takes a function, *priorityFn*, for computing these priority values.
The *setPriorities* function should be run after loading a complete set of data.
It will find each group of compounds which share the same
descriptor and call the given function, *priorityFn*,
with the compound_id numbers of the group. This function should
then assign priorities to each compound-descriptor pair, however
it wishes.
One built in priority function is *forestSizePriorities*. This simply
prefers compounds with fewer disconnected components over compounds with
more dissconnected components.
```{r eval=FALSE, tidy=FALSE}
setPriorities(conn,forestSizePriorities)
```
Searching
---------
Compounds can be searched for using the *findCompounds*
function. This function takes a connection object, a vector of feature
names used in the tests, and finally, a vector of tests that must all
pass for a compound to be included in the result set. Each test should
be a boolean expression. For example: 'c("MW <= 400","RINGS \> 3")'
would return all compounds with a molecular weight of 400 or less and
more than 3 rings, assuming these features exist in the database. The
syntax for each test is "\ \ \".
If you know SQL you can go beyond this basic syntax. These tests will
simply be concatenated together with "AND" in-between them and tacked on
the end of a WHERE clause of an SQL statement. So any SQL that will work
in that context is fine. The function will return a list of compound
ids, the actual compounds can be fetched with
*getCompounds*. If just the names are needed, the
*getCompoundNames* function can be used. Compounds can
also be fetched by name using the *findCompoundsByName*
function.
In this example we search for compounds with molecular weight less than
300.
```{r }
results = findCompounds(conn,"mw",c("mw < 300"))
message("found ",length(results))
```
If more than one test is given, only compounds which satisfy all tests are found. So if we wanted
to further restrict our search to compounds with 2 or more aromatic rings we could do:
```{r }
results = findCompounds(conn,c("mw","aromatic"),c("mw < 300","aromatic >= 2"))
message("found ",length(results))
```
Remember that any feature used in some test must be listed in the second argument.
String patterns can also be used. So if we wanted to match a substring of the molecular formula, say
to find compounds with 21 carbon atoms, we could do:
```{r }
results = findCompounds(conn,"formula",c("formula like '%C21%'"))
message("found ",length(results))
```
The "like" operator does a pattern match. There are two wildcard
operators that can be used with this operator. The "%" will match any stretch of characters while the "?"
will match any single character. So the above expression would match a formula like "C21H28N4O6".
Valid comparison operators are:
- <, <=, > , >=
- =, ==, !=, <>, IS, IS NOT, IN, LIKE
The boolean operators "AND" and "OR" can also be used to create more complex expressions within a single test.
If you just want to fetch every compound in the database you can use the *getAllCompoundIds* function:
```{r }
allIds = getAllCompoundIds(conn)
message("found ",length(allIds))
```
Using Search Results
-----------------------
Once you have a list of compound ids from the *findCompounds* function, you can either
fetch the compound names, or the whole set of compounds as an SDFset.
```{r }
#get the names of the compounds:
names = getCompoundNames(conn,results)
#if the name order is important set keepOrder=TRUE
#It will take a little longer though
names = getCompoundNames(conn,results,keepOrder=TRUE)
# get the whole set of compounds
compounds = getCompounds(conn,results)
#in order:
compounds = getCompounds(conn,results,keepOrder=TRUE)
#write results directly to a file:
compounds = getCompounds(conn,results,filename=file.path(tempdir(),"results.sdf"))
```
Using the *getCompoundFeatures* function, you can get a set of feature values
as a data frame:
```{r}
getCompoundFeatures(conn,results[1:5],c("mw","logp","formula"))
#write results directly to a CSV file (reduces memory usage):
getCompoundFeatures(conn,results[1:5],c("mw","logp","formula"),filename="features.csv")
#maintain input order in output:
print(results[1:5])
getCompoundFeatures(conn,results[1:5],c("mw","logp","formula"),keepOrder=TRUE)
```
Pre-Built Databases
--------------------
We have pre-built SQLite databases for the Drug Bank and DUD datasets. They can be found in
the ChemmineDrugs annotation package. Connections to these databases can be fetched from the
functions *DrugBank* and *DUD* to get the corresponding database. Any of the above functions can
then be used to query the database.
The DUD dataset was downloaded from [here](http://dude.docking.org/db/subsets/all/all.tar.gz). A description
can be found [here](http://dude.docking.org/).
The Drug Bank data set is version 4.1. It can be downloaded [here](http://www.drugbank.ca/system/downloads/current/structures/all.sdf.zip)
The following features are included:
- **aromatic**: Number of aromatic rings
- **cansmi**: Canonical SMILES sting
- **cansmins**:
- **formula**: Molecular formula
- **hba1**:
- **hba2**:
- **hbd**:
- **inchi**: INCHI string
- **logp**:
- **mr**:
- **mw**: Molecular weight
- **ncharges**:
- **nf**:
- **r2nh**:
- **r3n**:
- **rcch**:
- **rcho**:
- **rcn**:
- **rcooh**:
- **rcoor**:
- **rcor**:
- **rings**:
- **rnh2**:
- **roh**:
- **ropo3**:
- **ror**:
- **title**:
- **tpsa**:
The DUD database additionally includes:
- **target_name**: Name of the target
- **type**: either "active" or "decoy"
Working with SDF/SDFset Classes
===============================
Several methods are available to return the different data components of
*SDF/SDFset* containers in batches. The following
examples list the most important ones. To save space their content is
not printed in the manual.
```{r eval=FALSE, tidy=FALSE}
view(sdfset[1:4]) # Summary view of several molecules
length(sdfset) # Returns number of molecules
sdfset[[1]] # Returns single molecule from SDFset as SDF object
sdfset[[1]][[2]] # Returns atom block from first compound as matrix
sdfset[[1]][[2]][1:4,]
c(sdfset[1:4], sdfset[5:8]) # Concatenation of several SDFsets
```
The *grepSDFset* function allows string
matching/searching on the different data components in
*SDFset*. By default the function returns a SDF summary
of the matching entries. Alternatively, an index of the matches can be
returned with the setting *mode="index"*.
```{r eval=FALSE, tidy=FALSE}
grepSDFset("650001", sdfset, field="datablock", mode="subset") # To return index, set mode="index")
```
Utilities to maintain unique compound IDs:
```{r eval=FALSE, tidy=FALSE}
sdfid(sdfset[1:4]) # Retrieves CMP IDs from Molecule Name field in header block.
cid(sdfset[1:4]) # Retrieves CMP IDs from ID slot in SDFset.
unique_ids <- makeUnique(sdfid(sdfset)) # Creates unique IDs by appending a counter to duplicates.
cid(sdfset) <- unique_ids # Assigns uniquified IDs to ID slot
```
Subsetting by character, index and logical vectors:
```{r eval=FALSE, tidy=FALSE}
view(sdfset[c("650001", "650012")])
view(sdfset[4:1])
mylog <- cid(sdfset)
view(sdfset[mylog])
```
Accessing *SDF/SDFset* components: header, atom, bond and
data blocks:
```{r eval=FALSE, tidy=FALSE}
atomblock(sdf); sdf[[2]];
sdf[["atomblock"]] # All three methods return the same component
header(sdfset[1:4])
atomblock(sdfset[1:4])
bondblock(sdfset[1:4])
datablock(sdfset[1:4])
header(sdfset[[1]])
atomblock(sdfset[[1]])
bondblock(sdfset[[1]])
datablock(sdfset[[1]])
```
Replacement Methods:
```{r eval=FALSE, tidy=FALSE}
sdfset[[1]][[2]][1,1] <- 999
atomblock(sdfset)[1] <- atomblock(sdfset)[2]
datablock(sdfset)[1] <- datablock(sdfset)[2]
```
Assign matrix data to data block:
```{r eval=FALSE, tidy=FALSE}
datablock(sdfset) <- as.matrix(iris[1:100,])
view(sdfset[1:4])
```
Class coercions from *SDFstr* to *list*,
*SDF* and *SDFset*:
```{r eval=FALSE, tidy=FALSE}
as(sdfstr[1:2], "list") as(sdfstr[[1]], "SDF")
as(sdfstr[1:2], "SDFset")
```
Class coercions from *SDF* to *SDFstr*,
*SDFset*, list with SDF sub-components:
```{r eval=FALSE, tidy=FALSE}
sdfcomplist <- as(sdf, "list") sdfcomplist <-
as(sdfset[1:4], "list"); as(sdfcomplist[[1]], "SDF") sdflist <-
as(sdfset[1:4], "SDF"); as(sdflist, "SDFset") as(sdfset[[1]], "SDFstr")
as(sdfset[[1]], "SDFset")
```
Class coercions from *SDFset* to lists with components
consisting of SDF or sub-components:
```{r eval=FALSE, tidy=FALSE}
as(sdfset[1:4], "SDF") as(sdfset[1:4], "list") as(sdfset[1:4], "SDFstr")
```
Molecular Property Functions (Physicochemical Descriptors)
==========================================================
Several methods and functions are available to compute basic compound
descriptors, such as molecular formula (MF), molecular weight (MW), and
frequencies of atoms and functional groups. In many of these functions,
it is important to set *addH=TRUE* in order to
include/add hydrogens that are often not specified in an SD file.
```{r boxplot, eval=TRUE, tidy=FALSE}
propma <- atomcountMA(sdfset, addH=FALSE)
boxplot(propma, col="blue", main="Atom Frequency")
```
```{r eval=FALSE, tidy=FALSE}
boxplot(rowSums(propma), main="All Atom Frequency")
```
Data frame provided by library containing atom names, atom symbols,
standard atomic weights, group and period numbers:
```{r eval=TRUE, tidy=FALSE}
data(atomprop)
atomprop[1:4,]
```
Compute MW and formula:
```{r eval=TRUE, tidy=FALSE}
MW(sdfset[1:4], addH=FALSE)
MF(sdfset[1:4], addH=FALSE)
```
Enumerate functional groups:
```{r eval=TRUE, tidy=FALSE}
groups(sdfset[1:4], groups="fctgroup", type="countMA")
```
Combine MW, MF, charges, atom counts, functional group counts and ring
counts in one data frame:
```{r eval=TRUE, tidy=FALSE}
propma <- data.frame(MF=MF(sdfset, addH=FALSE), MW=MW(sdfset, addH=FALSE),
Ncharges=sapply(bonds(sdfset, type="charge"), length),
atomcountMA(sdfset, addH=FALSE),
groups(sdfset, type="countMA"),
rings(sdfset, upper=6, type="count", arom=TRUE))
propma[1:4,]
```
The following shows an example for assigning the values stored in a
matrix (*e.g.* property descriptors) to the data block components in an
*SDFset*. Each matrix row will be assigned to the
corresponding slot position in the *SDFset*.
```{r eval=FALSE, tidy=FALSE}
datablock(sdfset) <- propma # Works with all SDF components
datablock(sdfset)[1:4]
test <- apply(propma[1:4,], 1, function(x)
data.frame(col=colnames(propma), value=x))
```
The data blocks in SDFs contain often important annotation information
about compounds. The *datablock2ma* function returns this
information as matrix for all compounds stored in an
*SDFset* container. The *splitNumChar*
function can then be used to organize all numeric columns in a
*numeric matrix* and the character columns in a
*character matrix* as components of a
*list* object.
```{r eval=FALSE, tidy=FALSE}
datablocktag(sdfset, tag="PUBCHEM_NIST_INCHI")
datablocktag(sdfset,
tag="PUBCHEM_OPENEYE_CAN_SMILES")
```
Convert entire data block to matrix:
```{r eval=FALSE, tidy=FALSE}
blockmatrix <- datablock2ma(datablocklist=datablock(sdfset)) # Converts data block to matrix
numchar <- splitNumChar(blockmatrix=blockmatrix) # Splits matrix to numeric matrix and character matrix
numchar[[1]][1:4,]; numchar[[2]][1:4,]
# Splits matrix to numeric matrix and character matrix
```
Bond Matrices
=============
Bond matrices provide an efficient data structure for many basic
computations on small molecules. The function *conMA*
creates this data structure from *SDF* and
*SDFset* objects. The resulting bond matrix contains the
atom labels in the row/column titles and the bond types in the data
part. The labels are defined as follows: 0 is no connection, 1 is a
single bond, 2 is a double bond and 3 is a triple bond.
```{r contable, eval=FALSE, fig.keep='none', tidy=FALSE}
conMA(sdfset[1:2],
exclude=c("H")) # Create bond matrix for first two molecules in sdfset
conMA(sdfset[[1]], exclude=c("H")) # Return bond matrix for first molecule
plot(sdfset[1], atomnum = TRUE, noHbonds=FALSE , no_print_atoms = "", atomcex=0.8) # Plot its structure with atom numbering
rowSums(conMA(sdfset[[1]], exclude=c("H"))) # Return number of non-H bonds for each atom
```
Charges and Missing Hydrogens
=============================
The function *bonds* returns information about the number
of bonds, charges and missing hydrogens in *SDF* and
*SDFset* objects. It is used by many other functions
(*e.g.* *MW*, *MF*,
*atomcount*, *atomcuntMA* and
*plot*) to correct for missing hydrogens that are often
not specified in SD files.
```{r eval=TRUE, tidy=FALSE}
bonds(sdfset[[1]], type="bonds")[1:4,]
bonds(sdfset[1:2], type="charge")
bonds(sdfset[1:2], type="addNH")
```
Ring Perception and Aromaticity Assignment
==========================================
The function *rings* identifies all possible rings in one
or many molecules (here *sdfset[1]*) using the exhaustive
ring perception algorithm from `r citet(c("Hanser_1996"))` . In addition, the function can return
all smallest possible rings as well as aromaticity information.
The following example returns all possible rings in a
*list*. The argument *upper* allows to
specify an upper length limit for rings. Choosing smaller length limits
will reduce the search space resulting in shortened compute times. Note:
each ring is represented by a character vector of atom symbols that are
numbered by their position in the atom block of the corresponding
*SDF/SDFset* object.
```{r eval=TRUE, tidy=FALSE}
ringatoms <- rings(sdfset[1], upper=Inf, type="all", arom=FALSE, inner=FALSE)
```
For visual inspection, the corresponding compound structure can be
plotted with the ring bonds highlighted in color:
```{r eval=TRUE, tidy=FALSE}
atomindex <- as.numeric(gsub(".*_", "", unique(unlist(ringatoms))))
plot(sdfset[1], print=FALSE, colbonds=atomindex)
```
Alternatively, one can include the atom numbers in the plot:
```{r eval=TRUE, tidy=FALSE}
plot(sdfset[1], print=FALSE, atomnum=TRUE, no_print_atoms="H")
```
Aromaticity information of the rings can be returned in a logical vector
by setting *arom=TRUE*:
```{r eval=TRUE, tidy=FALSE}
rings(sdfset[1], upper=Inf, type="all", arom=TRUE, inner=FALSE)
```
Return rings with no more than 6 atoms that are also aromatic:
```{r eval=TRUE, tidy=FALSE}
rings(sdfset[1], upper=6, type="arom", arom=TRUE, inner=FALSE)
```
Count shortest possible rings and their aromaticity assignments by
setting *type=count* and *inner=TRUE*. The
inner (smallest possible) rings are identified by first computing all
possible rings and then selecting only the inner rings. For more
details, consult the help documentation with *?rings*.
```{r eval=TRUE, tidy=FALSE}
rings(sdfset[1:4], upper=Inf, type="count", arom=TRUE, inner=TRUE)
```
Rendering Chemical Structure Images
===================================
R Graphics Device
-----------------
A new plotting function for compound structures has been added to the
package recently. This function uses the native R graphics device for
generating compound depictions. At this point this function is still in
an experimental developmental stage but should become stable soon.
If you have *ChemmineOB* available you can use the *regenCoords*
option to have OpenBabel regenerate the coordinates for the compound.
This can sometimes produce better looking plots.
Plot compound Structures with R's graphics device:
```{r plotstruct2, eval=TRUE, tidy=FALSE}
data(sdfsample)
sdfset <- sdfsample
plot(sdfset[1:4], regenCoords=TRUE,print=FALSE) # 'print=TRUE' returns SDF summaries
```
Customized plots:
```{r eval=FALSE, tidy=FALSE}
plot(sdfset[1:4], griddim=c(2,2), print_cid=letters[1:4], print=FALSE,
noHbonds=FALSE)
```
In the following plot, the atom block position numbers in the SDF are
printed next to the atom symbols (*atomnum = TRUE*). For
more details, consult help documentation with
*?plotStruc* or *?plot*.
```{r plotstruct3, eval=TRUE, tidy=FALSE}
plot(sdfset["CMP1"], atomnum = TRUE, noHbonds=F , no_print_atoms = "",
atomcex=0.8, sub=paste("MW:", MW(sdfsample["CMP1"])), print=FALSE)
```
Substructure highlighting by atom numbers:
```{r plotstruct4, eval=TRUE, tidy=FALSE}
plot(sdfset[1], print=FALSE, colbonds=c(22,26,25,3,28,27,2,23,21,18,8,19,20,24))
```
Online with ChemMine Tools
--------------------------
Alternatively, one can visualize compound structures with a standard web
browser using the online ChemMine Tools service.
Plot structures using web service ChemMine Tools:
```{r eval=FALSE, tidy=FALSE}
sdf.visualize(sdfset[1:4])
```
Similarity Comparisons and Searching
====================================
Maximum Common Substructure (MCS) Searching
-------------------------------------------
The ChemmineR add-on package
[](http://www.bioconductor.org/packages/devel/bioc/html/fmcsR.html)
provides support for identifying maximum common substructures (MCSs) and
flexible MCSs among compounds. The algorithm can be used for pairwise
compound comparisons, structure similarity searching and clustering. The
manual describing this functionality is available
[](http://www.bioconductor.org/packages/devel/bioc/vignettes/fmcsR/inst/doc/fmcsR.pdf)
and the associated publication is available in
`r citet(c("Wang_2013"))`. The following gives a
short preview of some functionalities provided by the
*fmcsR* package.
```{r plotmcs, eval=TRUE, tidy=FALSE}
library(fmcsR)
data(fmcstest) # Loads test sdfset object
test <- fmcs(fmcstest[1], fmcstest[2], au=2, bu=1) # Searches for MCS with mismatches
plotMCS(test) # Plots both query compounds with MCS in color
```
AP/APset Classes for Storing Atom Pair Descriptors
--------------------------------------------------
The function *sdf2ap* computes atom pair descriptors for
one or many compounds
`r citep(c("Carhart_1985","Chen_2002"))` . It returns a searchable atom pair database
stored in a container of class *APset*, which can be used
for structural similarity searching and clustering. As similarity
measure, the Tanimoto coefficient or related coefficients can be used.
An *APset* object consists of one or many
*AP* entries each storing the atom pairs of a single
compound. Note: the deprecated *cmp.parse* function is
still available which also generates atom pair descriptor databases, but
directly from an SD file. Since the latter function is less flexible it
may be discontinued in the future.
Generate atom pair descriptor database for searching:
```{r eval=TRUE, tidy=FALSE}
ap <- sdf2ap(sdfset[[1]]) # For single compound
ap
```
```{r eval=FALSE, tidy=FALSE}
apset <- sdf2ap(sdfset)
# For many compounds.
```
```{r eval=TRUE, tidy=FALSE}
view(apset[1:4])
```
Return main components of APset objects:
```{r eval=FALSE, tidy=FALSE}
cid(apset[1:4]) # Compound IDs
ap(apset[1:4]) # Atom pair
descriptors
db.explain(apset[1]) # Return atom pairs in human readable format
```
Coerce APset to other objects:
```{r eval=FALSE, tidy=FALSE}
apset2descdb(apset) # Returns old list-style AP database
tmp <- as(apset, "list") # Returns list
as(tmp, "APset") # Converts list back to APset
```
Large SDF and Atom Pair Databases
---------------------------------
When working with large data sets it is often desirable to save the
*SDFset* and *APset* containers as binary
R objects to files for later use. This way they can be loaded very
quickly into a new R session without recreating them every time from
scratch.
Save and load of *SDFset* and *APset*
containers:
```{r eval=FALSE, tidy=FALSE}
save(sdfset, file = "sdfset.rda", compress = TRUE)
load("sdfset.rda") save(apset, file = "apset.rda", compress = TRUE)
load("apset.rda")
```
Pairwise Compound Comparisons with Atom Pairs
---------------------------------------------
The *cmp.similarity* function computes the atom pair
similarity between two compounds using the Tanimoto coefficient as
similarity measure. The coefficient is defined as *c/(a+b+c)*, which
is the proportion of the atom pairs shared among two compounds divided
by their union. The variable *c* is the number of atom pairs common in
both compounds, while *a* and *b* are the numbers of their unique
atom pairs.
```{r eval=TRUE, tidy=FALSE}
cmp.similarity(apset[1],
apset[2])
cmp.similarity(apset[1], apset[1])
```
Similarity Searching with Atom Pairs
------------------------------------
The *cmp.search* function searches an atom pair database
for compounds that are similar to a query compound. The following
example returns a data frame where the rows are sorted by the Tanimoto
similarity score (best to worst). The first column contains the indices
of the matching compounds in the database. The argument cutoff can be a
similarity cutoff, meaning only compounds with a similarity value larger
than this cutoff will be returned; or it can be an integer value
restricting how many compounds will be returned. When supplying a cutoff
of 0, the function will return the similarity values for every compound
in the database.
```{r eval=TRUE, tidy=FALSE}
cmp.search(apset,
apset["650065"], type=3, cutoff = 0.3, quiet=TRUE)
```
Alternatively, the
function can return the matches in form of an index or a named vector if
the *type* argument is set to *1* or
*2*, respectively.
```{r eval=TRUE, tidy=FALSE}
cmp.search(apset, apset["650065"], type=1, cutoff = 0.3, quiet=TRUE)
cmp.search(apset, apset["650065"], type=2, cutoff = 0.3, quiet=TRUE)
```
FP/FPset Classes for Storing Fingerprints
-----------------------------------------
The *FPset* class stores fingerprints of small molecules
in a matrix-like representation where every molecule is encoded as a
fingerprint of the same type and length. The *FPset*
container acts as a searchable database that contains the fingerprints
of many molecules. The *FP* container holds only one
fingerprint. Several constructor and coerce methods are provided to
populate *FP/FPset* containers with fingerprints, while
supporting any type and length of fingerprints. For instance, the
function *desc2fp* generates fingerprints from an atom
pair database stored in an *APset*, and
*as(matrix, "FPset")* and *as(character,
"FPset")* construct an *FPset* database from
objects where the fingerprints are represented as
*matrix* or *character* objects,
respectively.
Show slots of *FPset* class:
```{r eval=TRUE, tidy=FALSE}
showClass("FPset")
```
Instance of *FPset* class:
```{r eval=TRUE, tidy=FALSE}
data(apset)
fpset <- desc2fp(apset)
view(fpset[1:2])
```
*FPset* class usage:
```{r eval=TRUE, tidy=FALSE}
fpset[1:4] # behaves like a list
fpset[[1]] # returns FP object
length(fpset) # number of compounds ENDCOMMENT
cid(fpset) # returns compound ids
fpset[10] <- 0 # replacement of 10th fingerprint to all zeros
cid(fpset) <- 1:length(fpset) # replaces compound ids
c(fpset[1:4], fpset[11:14]) # concatenation of several FPset objects
```
Construct *FPset* class form *matrix*:
```{r eval=TRUE, tidy=FALSE}
fpma <- as.matrix(fpset) # coerces FPset to matrix
as(fpma, "FPset")
```
Construct *FPset* class form *character
vector*:
```{r eval=TRUE, tidy=FALSE}
fpchar <- as.character(fpset) # coerces FPset to character strings
as(fpchar, "FPset") # construction of FPset class from character vector
```
Compound similarity searching with *FPset*:
```{r eval=TRUE, tidy=FALSE}
fpSim(fpset[1], fpset, method="Tanimoto", cutoff=0.4, top=4)
```
Folding fingerprints:
```{r eval=TRUE,tidy=FALSE}
fold(fpset) # fold each FP once
fold(fpset, count=2) #fold each FP twice
fold(fpset, bits=128) #fold each FP down to 128 bits
fold(fpset[[1]]) # fold an individual FP
fptype(fpset) # get type of FPs
numBits(fpset) # get the number of bits of each FP
foldCount(fold(fpset)) # the number of times an FP or FPset has been folded
```
Atom Pair Fingerprints
----------------------
Atom pairs can be converted into binary atom pair fingerprints of fixed
length. Computations on this compact data structure are more time and
memory efficient than on their relatively complex atom pair
counterparts. The function *desc2fp* generates
fingerprints from descriptor vectors of variable length such as atom
pairs stored in *APset* or *list*
containers. The obtained fingerprints can be used for structure
similarity comparisons, searching and clustering.
Create atom pair sample data set:
```{r eval=FALSE, tidy=FALSE}
data(sdfsample)
sdfset <- sdfsample[1:10]
apset <- sdf2ap(sdfset)
```
Compute atom pair fingerprint database using internal atom pair
selection containing the 4096 most common atom pairs identified in
DrugBank's compound collection. For details see *?apfp*.
The following example uses from this set the 1024 most frequent atom
pairs:
```{r eval=FALSE, tidy=FALSE}
fpset <- desc2fp(apset, descnames=1024, type="FPset")
```
Alternatively, one can provide any custom atom pair selection. Here, the
1024 most common ones in *apset*:
```{r eval=FALSE, tidy=FALSE}
fpset1024 <- names(rev(sort(table(unlist(as(apset, "list")))))[1:1024])
fpset <- desc2fp(apset, descnames=fpset1024, type="FPset")
```
A more compact way of storing fingerprints is as character values:
```{r eval=FALSE, tidy=FALSE}
fpchar <- desc2fp(x=apset,
descnames=1024, type="character") fpchar <- as.character(fpset)
```
Converting a fingerprint database to a matrix and vice versa:
```{r eval=FALSE, tidy=FALSE}
fpma <- as.matrix(fpset)
fpset <- as(fpma, "FPset")
```
Similarity searching and returning Tanimoto similarity coefficients:
```{r eval=FALSE, tidy=FALSE}
fpSim(fpset[1], fpset, method="Tanimoto")
```
Under *method* one can choose from several predefined
similarity measures including *Tanimoto* (default),
*Euclidean*, *Tversky* or
*Dice*. Alternatively, one can pass on custom similarity
functions.
```{r eval=FALSE, tidy=FALSE}
fpSim(fpset[1], fpset, method="Tversky", cutoff=0.4, top=4, alpha=0.5, beta=1)
```
Example for using a custom similarity function:
```{r eval=FALSE, tidy=FALSE}
myfct <- function(a, b, c, d) c/(a+b+c+d)
fpSim(fpset[1], fpset, method=myfct)
```
Clustering example:
```{r eval=FALSE, tidy=FALSE}
simMAap <- sapply(cid(apfpset), function(x) fpSim(x=apfpset[x], apfpset, sorted=FALSE))
hc <- hclust(as.dist(1-simMAap), method="single")
plot(as.dendrogram(hc), edgePar=list(col=4, lwd=2), horiz=TRUE)
```
Fingerprint E-values
---------------------
The *fpSim* function can also return Z-scores, E-values, and p-values
if given a set of score distribution parameters. These parameters can
be computed over an *fpSet* with the *genParameters* function.
```{r eval=TRUE, tidy=FALSE}
params = genParameters(fpset)
```
This function will compute all pairwise distances between the given
fingerprints and then fit a Beta distribution to the resulting
Tanimoto scores, conditioned on the number of set bits in each
fingerprint. For large data sets where you would not want to compute
all pairwise distances, you can set what fraction to sample with the
*sampleFraction* argument. This step only needs to be done once for
each database of *fpSet* objects. Alternatively, if you have a large
database of fingerprints, or you believe that the parameters computed
on a single database are more generally applicable, you can use the
resulting parameters for other databases as well.
Once you have a set of parameters, you can pass them to *fpSim* with
the *parameters* argument.
```{r eval=TRUE, tidy=FALSE}
fpSim(fpset[[1]],fpset,top=10,parameters=params)
```
This will then return a data frame with the similarity, Z-score,
E-value, and p-value. You can change which value will be used as a
cutoff and to sort by by setting the argument *scoreType* to one of
these scores. In this way you could set an E-value cutoff of 0.04 for
example.
```{r eval=TRUE, tidy=FALSE}
fpSim(fpset[[1]],fpset,cutoff=0.04,scoreType="evalue",parameters=params)
```
Pairwise Compound Comparisons with PubChem Fingerprints
-------------------------------------------------------
The *fpSim* function computes the similarity coefficients
(*e.g.* Tanimoto) for pairwise comparisons of binary fingerprints. For
this data type, *c* is the number of "on-bits" common in both
compounds, and *a* and *b* are the numbers of their unique
"on-bits". Currently, the PubChem fingerprints need to be provided (here
PubChem's SD files) and cannot be computed from scratch in
*ChemmineR*. The PubChem fingerprint specifications can
be loaded with *data(pubchemFPencoding)*.
Convert base 64 encoded PubChem fingerprints to
*character* vector, *matrix* or
*FPset* object:
```{r eval=TRUE, tidy=FALSE}
cid(sdfset) <- sdfid(sdfset)
fpset <- fp2bit(sdfset, type=1)
fpset <- fp2bit(sdfset, type=2)
fpset <- fp2bit(sdfset, type=3)
fpset
```
Pairwise compound structure comparisons:
```{r eval=TRUE, tidy=FALSE}
fpSim(fpset[1], fpset[2])
```
Similarity Searching with PubChem Fingerprints
----------------------------------------------
Similarly, the *fpSim* function provides search
functionality for PubChem fingerprints:
```{r eval=TRUE, tidy=FALSE}
fpSim(fpset["650065"], fpset, method="Tanimoto", cutoff=0.6, top=6)
```
Visualize Similarity Search Results
-----------------------------------
The *cmp.search* function allows to visualize the
chemical structures for the search results. Similar but more flexible
chemical structure rendering functions are *plot* and
*sdf.visualize* described above. By setting the visualize
argument in *cmp.search* to *TRUE*, the
matching compounds and their scores can be visualized with a standard
web browser. Depending on the *visualize.browse*
argument, an URL will be printed or a webpage will be opened showing the
structures of the matching compounds.
View similarity search results in R's graphics device:
```{r search_result, eval=TRUE, tidy=FALSE}
cid(sdfset) <-
cid(apset) # Assure compound name consistency among objects.
plot(sdfset[names(cmp.search(apset, apset["650065"], type=2, cutoff=4, quiet=TRUE))], print=FALSE)
```
View results online with Chemmine Tools:
```{r eval=FALSE, tidy=FALSE}
similarities <- cmp.search(apset, apset[1], type=3, cutoff = 10)
sdf.visualize(sdfset[similarities[,1]])
```
Clustering
==========
Clustering Identical or Very Similar Compounds
----------------------------------------------
Often it is of interest to identify very similar or identical compounds
in a compound set. The *cmp.duplicated* function can be
used to quickly identify very similar compounds in atom pair sets, which
will be frequently, but not necessarily, identical compounds.
Identify compounds with identical AP sets:
```{r eval=TRUE, tidy=FALSE}
cmp.duplicated(apset, type=1)[1:4] # Returns AP duplicates as logical vector
cmp.duplicated(apset, type=2)[1:4,] # Returns AP duplicates as data frame
```
Plot the structure of two pairs of duplicates:
```{r duplicates, eval=TRUE, tidy=FALSE}
plot(sdfset[c("650059","650060", "650065", "650066")], print=FALSE)
```
Remove AP duplicates from SDFset and APset objects:
```{r eval=TRUE, tidy=FALSE}
apdups <- cmp.duplicated(apset, type=1)
sdfset[which(!apdups)]; apset[which(!apdups)]
```
Alternatively, one can identify duplicates via other descriptor types if
they are provided in the data block of an imported SD file. For
instance, one can use here fingerprints, InChI, SMILES or other
molecular representations. The following examples show how to enumerate
by identical InChI strings, SMILES strings and molecular formula,
respectively.
```{r eval=TRUE, tidy=FALSE}
count <- table(datablocktag(sdfset,
tag="PUBCHEM_NIST_INCHI"))
count <- table(datablocktag(sdfset, tag="PUBCHEM_OPENEYE_CAN_SMILES"))
count <- table(datablocktag(sdfset, tag="PUBCHEM_MOLECULAR_FORMULA"))
count[1:4]
```
Binning Clustering
------------------
Compound libraries can be clustered into discrete similarity groups with
the binning clustering function *cmp.cluster*. The
function accepts as input an atom pair (*APset*) or a
fingerprint (*FPset*) descriptor database as well as a
similarity threshold. The binning clustering result is returned in form
of a data frame. Single linkage is used for cluster joining. The
function calculates the required compound-to-compound distance
information on the fly, while a memory-intensive distance matrix is only
created upon user request via the *save.distances*
argument (see below).
Because an optimum similarity threshold is often not known, the
*cmp.cluster* function can calculate cluster results for
multiple cutoffs in one step with almost the same speed as for a single
cutoff. This can be achieved by providing several cutoffs under the
cutoff argument. The clustering results for the different cutoffs will
be stored in one data frame.
One may force the *cmp.cluster* function to calculate and
store the distance matrix by supplying a file name to the
*save.distances* argument. The generated distance matrix
can be loaded and passed on to many other clustering methods available
in R, such as the hierarchical clustering function
*hclust* (see below).
If a distance matrix is available, it may also be supplied to
*cmp.cluster* via the *use.distances*
argument. This is useful when one has a pre-computed distance matrix
either from a previous call to *cmp.cluster* or from
other distance calculation subroutines.
Single-linkage binning clustering with one or multiple cutoffs:
```{r eval=TRUE, tidy=FALSE}
clusters <- cmp.cluster(db=apset, cutoff = c(0.7, 0.8, 0.9), quiet = TRUE)
clusters[1:12,]
```
Clustering of *FPset* objects with multiple cutoffs. This
method allows to call various similarity methods provided by the
*fpSim* function. For details consult
*?fpSim*.
```{r eval=TRUE, tidy=FALSE}
fpset <- desc2fp(apset)
clusters2 <- cmp.cluster(fpset, cutoff=c(0.5, 0.7, 0.9), method="Tanimoto", quiet=TRUE)
clusters2[1:12,]
```
Sames as above, but using Tversky similarity measure:
```{r eval=TRUE, tidy=FALSE}
clusters3 <- cmp.cluster(fpset, cutoff=c(0.5, 0.7, 0.9),
method="Tversky", alpha=0.3, beta=0.7, quiet=TRUE)
```
Return cluster size distributions for each cutoff:
```{r eval=TRUE, tidy=FALSE}
cluster.sizestat(clusters, cluster.result=1)
cluster.sizestat(clusters, cluster.result=2)
cluster.sizestat(clusters, cluster.result=3)
```
Enforce calculation of distance matrix:
```{r eval=FALSE, tidy=FALSE}
clusters <- cmp.cluster(db=apset, cutoff = c(0.65, 0.5, 0.3),
save.distances="distmat.rda") # Saves distance matrix to file "distmat.rda" in current working directory.
load("distmat.rda") # Loads distance matrix.
```
Jarvis-Patrick Clustering
-------------------------
The Jarvis-Patrick clustering algorithm is widely used in
cheminformatics `r citep("greycite13371")`. It requires a nearest neighbor table, which consists
of *j* nearest neighbors for each item (*e.g.* compound).
The nearest neighbor table is then used to join items into clusters when
they meet the following requirements: (a) they are contained in each
other's neighbor list and (b) they share at least *k*
nearest neighbors. The values for *j* and
*k* are user-defined parameters. The
*jarvisPatrick* function implemented in
*ChemmineR* takes a nearest neighbor table generated by
*nearestNeighbors*, which works for
*APset* and *FPset* objects. This function
takes either the standard Jarvis-Patrick *j* parameter
(as the *numNbrs* parameter), or else a
*cutoff* value, which is an extension to the basic
algorithm that we have added. Given a cutoff value, the nearest neighbor
table returned contains every neighbor with a similarity greater than
the cutoff value, for each item. This allows one to generate tighter
clusters and to minimize certain limitations of this method, such as
false joins of completely unrelated items when operating on small data
sets. The *trimNeighbors* function can also be used to
take an existing nearest neighbor table and remove all neighbors whose
similarity value is below a given cutoff value. This allows one to
compute a very relaxed nearest neighbor table initially, and then
quickly try different refinements later.
In case an existing nearest neighbor matrix needs to be used, the
*fromNNMatrix* function can be used to transform it into
the list structure that *jarvisPatrick* requires. The
input matrix must have a row for each compound, and each row should be
the index values of the neighbors of compound represented by that row.
The names of each compound can also be given through the
*names* argument. If not given, it will attempt to use
the *rownames* of the given matrix.
The *jarvisPatrick* function also allows one to relax
some of the requirements of the algorithm through the
*mode* parameter. When set to "a1a2b", then all
requirements are used. If set to "a1b", then (a) is relaxed to a
unidirectional requirement. Lastly, if *mode* is set to
"b", then only requirement (b) is used, which means that all pairs of
items will be checked to see if (b) is satisfied between them. The size
of the clusters generated by the different methods increases in this
order: "a1a2b" < "a1b" < "b". The run time of method "a1a2b" follows a
close to linear relationship, while it is nearly quadratic for the much
more exhaustive method "b". Only methods "a1a2b" and "a1b" are suitable
for clustering very large data sets (e.g. \>50,000 items) in a
reasonable amount of time.
An additional extension to the algorithm is the ability to set the
linkage mode. The *linkage* parameter can be one of
"single", "average", or "complete", for single linkage, average linkage
and complete linkage merge requirements, respectively. In the context of
Jarvis-Patrick, average linkage means that at least half of the pairs
between the clusters under consideration must meet requirement (b).
Similarly, for complete linkage, all pairs must requirement (b). Single
linkage is the normal case for Jarvis-Patrick and just means that at
least one pair must meet requirement (b).
The output is a cluster *vector* with the item labels in
the name slot and the cluster IDs in the data slot. There is a utility
function called *byCluster*, which takes out cluster
vector output by *jarvisPatrick* and transforms it into a
list of vectors. Each slot of the list is named with a cluster id and
the vector contains the cluster members. By default the function
excludes singletons from the output, but they can be included by setting
*excludeSingletons*=FALSE.
Load/create sample *APset* and *FPset*:
```{r eval=TRUE, tidy=FALSE}
data(apset)
fpset <- desc2fp(apset)
```
Standard Jarvis-Patrick clustering on *APset* and
*FPset* objects:
```{r eval=TRUE, tidy=FALSE}
jarvisPatrick(nearestNeighbors(apset,numNbrs=6), k=5, mode="a1a2b")
#Using "APset"
jarvisPatrick(nearestNeighbors(fpset,numNbrs=6), k=5, mode="a1a2b")
#Using "FPset"
```
The following example runs Jarvis-Patrick clustering with a minimum
similarity *cutoff* value (here Tanimoto coefficient). In
addition, it uses the much more exhaustive *"b"* method
that generates larger cluster sizes, but significantly increased the run
time. For more details, consult the corresponding help file with
*?jarvisPatrick*.
```{r eval=TRUE, tidy=FALSE}
cl<-jarvisPatrick(nearestNeighbors(fpset,cutoff=0.6,
method="Tanimoto"), k=2 ,mode="b")
byCluster(cl)
```
Output nearest neighbor table (*matrix*):
```{r eval=TRUE, tidy=FALSE}
nnm <- nearestNeighbors(fpset,numNbrs=6)
nnm$names[1:4]
nnm$ids[1:4,]
nnm$similarities[1:4,]
```
Trim nearest neighbor table:
```{r eval=TRUE, tidy=FALSE}
nnm <- trimNeighbors(nnm,cutoff=0.4)
nnm$similarities[1:4,]
```
Perform clustering on precomputed nearest neighbor table:
```{r eval=TRUE, tidy=FALSE}
jarvisPatrick(nnm, k=5,mode="b")
```
Using a user defined nearest neighbor matrix:
```{r eval=TRUE, tidy=FALSE}
nn <- matrix(c(1,2,2,1),2,2,dimnames=list(c('one','two')))
nn
byCluster(jarvisPatrick(fromNNMatrix(nn),k=1))
```
Multi-Dimensional Scaling (MDS)
-------------------------------
To visualize and compare clustering results, the
*cluster.visualize* function can be used. The function
performs Multi-Dimensional Scaling (MDS) and visualizes the results in
form of a scatter plot. It requires as input an *APset*,
a clustering result from *cmp.cluster*, and a cutoff for
the minimum cluster size to consider in the plot. To help determining a
proper cutoff size, the *cluster.sizestat* function is
provided to generate cluster size statistics.
MDS clustering and scatter plot:
```{r eval=FALSE, tidy=FALSE}
cluster.visualize(apset, clusters, size.cutoff=2, quiet = TRUE) # Color codes clusters with at least two members.
cluster.visualize(apset, clusters, quiet = TRUE) # Plots all items.
```
Create a 3D scatter plot of MDS result:
```{r mds_scatter, eval=TRUE, tidy=FALSE}
library(scatterplot3d)
coord <- cluster.visualize(apset, clusters, size.cutoff=1, dimensions=3, quiet=TRUE)
scatterplot3d(coord)
```
Interactive 3D scatter plot with Open GL (graphics not evaluated here):
```{r eval=FALSE, tidy=FALSE}
library(rgl) rgl.open(); offset <- 50;
par3d(windowRect=c(offset, offset, 640+offset, 640+offset))
rm(offset)
rgl.clear()
rgl.viewpoint(theta=45, phi=30, fov=60, zoom=1)
spheres3d(coord[,1], coord[,2], coord[,3], radius=0.03, color=coord[,4], alpha=1, shininess=20)
aspect3d(1, 1, 1)
axes3d(col='black')
title3d("", "", "", "", "", col='black')
bg3d("white") # To save a snapshot of the graph, one can use the command rgl.snapshot("test.png").
```
Clustering with Other Algorithms
--------------------------------
*ChemmineR* allows the user to take advantage of the wide
spectrum of clustering utilities available in R. An example on how to
perform hierarchical clustering with the hclust function is given
below.
Create atom pair distance matrix:
```{r ap_dist_matrix, eval=TRUE, tidy=FALSE}
dummy <- cmp.cluster(db=apset, cutoff=0, save.distances="distmat.rda", quiet=TRUE)
load("distmat.rda")
```
Hierarchical clustering with *hclust*:
```{r hclust, eval=TRUE, tidy=FALSE}
hc <- hclust(as.dist(distmat), method="single")
hc[["labels"]] <- cid(apset) # Assign correct item labels
plot(as.dendrogram(hc), edgePar=list(col=4, lwd=2), horiz=T)
```
Instead of atom pairs one can use PubChem's fingerprints for clustering:
```{r fp_hclust, eval=FALSE, tidy=FALSE}
simMA <- sapply(cid(fpset), function(x) fpSim(fpset[x], fpset, sorted=FALSE))
hc <- hclust(as.dist(1-simMA), method="single")
plot(as.dendrogram(hc), edgePar=list(col=4, lwd=2), horiz=TRUE)
```
Plot dendrogram with heatmap (here similarity matrix):
```{r heatmap, eval=TRUE, tidy=FALSE}
library(gplots)
heatmap.2(1-distmat, Rowv=as.dendrogram(hc), Colv=as.dendrogram(hc),
col=colorpanel(40, "darkblue", "yellow", "white"),
density.info="none", trace="none")
```
Searching PubChem
=================
Get Compounds from PubChem by Id
--------------------------------
The function *getIds* accepts one or more numeric PubChem
compound ids and downloads the corresponding compounds from PubChem
Power User Gateway (PUG) returning results in an *SDFset*
container. The ChemMine Tools web service is used as an intermediate, to
translate queries from plain HTTP POST to a PUG SOAP query.
Fetch 2 compounds from PubChem:
```{r getIds, eval=FALSE, tidy=FALSE}
compounds <- getIds(c(111,123))
compounds
```
Search a SMILES Query in PubChem
--------------------------------
The function *searchString* accepts one SMILES string
(Simplified Molecular Input Line Entry Specification) and performs a
\>0.95 similarity PubChem fingerprint search, returning the hits in an
*SDFset* container. The ChemMine Tools web service is
used as an intermediate, to translate queries from plain HTTP POST to a
PubChem Power User Gateway (PUG) query.
Search a SMILES string on PubChem:
```{r searchString, eval=FALSE, tidy=FALSE}
compounds <- searchString("CC(=O)OC1=CC=CC=C1C(=O)O") compounds
```
Search an SDF Query in PubChem
------------------------------
The function *searchSim* performs a PubChem similarity
search just like *searchString*, but accepts a query in
an *SDFset* container. If the query contains more than
one compound, only the first is searched.
Search an *SDFset* container on PubChem:
```{r searchSim, eval=FALSE, tidy=FALSE}
data(sdfsample);
sdfset <- sdfsample[1]
compounds <- searchSim(sdfset)
compounds
```
ChemMine Tools R Interface
==========================
ChemMine Web Tools is an online service for analyzing and clustering small molecules. It provides numerous cheminformatics tools which can be used directly on the website, or called remotely from within R. When called within R jobs are sent remotely to a queue on a compute cluster at UC Riverside, which is a free service offered to *ChemmineR* users.
The website is free and open to all users and is available at . When new tools are added to the service, they automatically become availiable within *ChemmineR* without updating your local R package.
List all available tools:
```{r listCMTools, eval=FALSE, tidy=FALSE}
listCMTools()
```
Show options and description for a tool. This also provides an example function call which can be copied
verbatim, and changed as necessary:
```{r toolDetailsCMT, eval=TRUE, tidy=FALSE}
toolDetails("Fingerprint Search")
```
Launch a Job
------------------------------
When a job is launched it returns a job token which refers to the running job on the UC Riverside cluster. You can check the status of a job or obtain the results as follows. If *result* is called on a job that is still running, it will loop internally until the job is completed, and then return the result.
Launch the tool *pubchemID2SDF* to obtain the structure for PubChem cid 2244:
```{r launchCMTool, eval=FALSE, tidy=FALSE}
job1 <- launchCMTool("pubchemID2SDF", 2244)
status(job1)
result1 <- result(job1)
```
Use the previous result to search PubChem for similar compounds:
```{r fingerprintSearchCMT, eval=FALSE, tidy=FALSE}
job2 <- launchCMTool('Fingerprint Search', result1, 'Similarity Cutoff'=0.95, 'Max Compounds Returned'=200)
result2 <- result(job2)
job3 <- launchCMTool("pubchemID2SDF", result2)
result3 <- result(job3)
```
Compute OpenBabel descriptors for these search results:
```{r obDescriptorsCMT, eval=FALSE, tidy=FALSE}
job4 <- launchCMTool("OpenBabel Descriptors", result3)
result4 <- result(job4)
result4
```
View Job Result in Browser
------------------------------
The function *browseJob* launches a web browser to view the results of a job online, just as if they
had been run from the ChemMine Tools website itself. If you also want the result data within R, you must first call
the *result* object from within R before calling *browseJob*. Once *browseJob* has been called on a job token,
the results are no longer accessible within R.
If you have an account on ChemMine Tools and would like to save the web results from your job, you must first login to your account within the default web browser on your system before you launch *browseJob*. The job will then be assigned automatically to the currently logged in account.
View OpenBabel descriptors online:
```{r obDescriptorsWWW, eval=FALSE, tidy=FALSE}
browseJob(job4)
```
Perform binning clustering and visualize result online:
```{r binningClusterWWW, eval=FALSE, tidy=FALSE}
job5 <- launchCMTool("Binning Clustering", result3, 'Similarity Cutoff'=0.9)
browseJob(job5)
```
Version Information
===================
```{r sessionInfo, results='asis'}
sessionInfo()
```
Funding
=======
This software was developed with funding from the National Science
Foundation:
[](http://www.nsf.gov/awardsearch/showAward.do?AwardNumber=0957099),
2010-0520325 and IGERT-0504249.
References
===========
```{r biblio, echo=FALSE, results='asis'}
#hack to get rid of latex code returned by certain doi sites (e.g. crosscite.org)
fix = function(t) gsub("\\$\\\\(l|r)brace\\$", "", t)
suppressWarnings({
sink(tempfile())
b<<-bibliography()
sink(NULL)
})
for(i in 1:length(b))
b[i]$title = fix(b[i]$title)
# write.bibtex(file="references.bib")
b
```