Source code

001/*
002 *                    BioJava development code
003 *
004 * This code may be freely distributed and modified under the
005 * terms of the GNU Lesser General Public Licence.  This should
006 * be distributed with the code.  If you do not have a copy,
007 * see:
008 *
009 *      http://www.gnu.org/copyleft/lesser.html
010 *
011 * Copyright for this code is held jointly by the individual
012 * authors.  These should be listed in @author doc comments.
013 *
014 * For more information on the BioJava project and its aims,
015 * or to join the biojava-l mailing list, visit the home page
016 * at:
017 *
018 *      http://www.biojava.org/
019 *
020 * Created on June 24, 2010
021 * Author: Mark Chapman
022 */
023
024package org.biojava.nbio.alignment.template;
025
026import org.biojava.nbio.core.alignment.template.ProfilePair;
027import org.biojava.nbio.core.alignment.template.Profile;
028import org.biojava.nbio.core.alignment.template.SubstitutionMatrix;
029import org.biojava.nbio.alignment.template.GapPenalty.Type;
030import org.biojava.nbio.core.sequence.template.Compound;
031import org.biojava.nbio.core.sequence.template.CompoundSet;
032import org.biojava.nbio.core.sequence.template.Sequence;
033import org.slf4j.Logger;
034import org.slf4j.LoggerFactory;
035
036import java.util.ArrayList;
037import java.util.List;
038import java.util.concurrent.ExecutionException;
039import java.util.concurrent.Future;
040
041/**
042 * Implements common code for an {@link Aligner} for a pair of {@link Profile}s.
043 *
044 * @author Mark Chapman
045 * @param <S> each {@link Sequence} in the pair of alignment {@link Profile}s is of type S
046 * @param <C> each element of an {@link AlignedSequence} is a {@link Compound} of type C
047 */
048public abstract class AbstractProfileProfileAligner<S extends Sequence<C>, C extends Compound>
049                extends AbstractMatrixAligner<S, C> implements ProfileProfileAligner<S, C> {
050
051        private final static Logger logger = LoggerFactory.getLogger(AbstractProfileProfileAligner.class);
052
053        // additional input fields
054        private Profile<S, C> query, target;
055
056        // concurrent execution fields
057        private Future<ProfilePair<S, C>> queryFuture, targetFuture;
058
059        // cached fields
060        private List<C> cslist;
061        private float[][] qfrac, tfrac;
062
063        // additional output field
064        protected ProfilePair<S, C> pair;
065
066        /**
067         * Before running a profile-profile alignment, data must be sent in via calls to
068         * {@link #setQuery(Profile)}, {@link #setTarget(Profile)}, {@link #setGapPenalty(GapPenalty)}, and
069         * {@link #setSubstitutionMatrix(SubstitutionMatrix)}.
070         */
071        protected AbstractProfileProfileAligner() {
072        }
073
074        /**
075         * Prepares for a profile-profile alignment.
076         *
077         * @param query the first {@link Profile} of the pair to align
078         * @param target the second {@link Profile} of the pair to align
079         * @param gapPenalty the gap penalties used during alignment
080         * @param subMatrix the set of substitution scores used during alignment
081         */
082        protected AbstractProfileProfileAligner(Profile<S, C> query, Profile<S, C> target, GapPenalty gapPenalty,
083                        SubstitutionMatrix<C> subMatrix) {
084                super(gapPenalty, subMatrix);
085                this.query = query;
086                this.target = target;
087                reset();
088        }
089
090        /**
091         * Prepares for a profile-profile alignment run concurrently.
092         *
093         * @param query the first {@link Profile} of the pair to align, still to be calculated
094         * @param target the second {@link Profile} of the pair to align, still to be calculated
095         * @param gapPenalty the gap penalties used during alignment
096         * @param subMatrix the set of substitution scores used during alignment
097         */
098        protected AbstractProfileProfileAligner(Future<ProfilePair<S, C>> query, Future<ProfilePair<S, C>> target,
099                        GapPenalty gapPenalty, SubstitutionMatrix<C> subMatrix) {
100                super(gapPenalty, subMatrix);
101                queryFuture = query;
102                targetFuture = target;
103                reset();
104        }
105
106        /**
107         * Prepares for a profile-profile alignment run concurrently.
108         *
109         * @param query the first {@link Profile} of the pair to align
110         * @param target the second {@link Profile} of the pair to align, still to be calculated
111         * @param gapPenalty the gap penalties used during alignment
112         * @param subMatrix the set of substitution scores used during alignment
113         */
114        protected AbstractProfileProfileAligner(Profile<S, C> query, Future<ProfilePair<S, C>> target,
115                        GapPenalty gapPenalty, SubstitutionMatrix<C> subMatrix) {
116                super(gapPenalty, subMatrix);
117                this.query = query;
118                targetFuture = target;
119                reset();
120        }
121
122        /**
123         * Prepares for a profile-profile alignment run concurrently.
124         *
125         * @param query the first {@link Profile} of the pair to align, still to be calculated
126         * @param target the second {@link Profile} of the pair to align
127         * @param gapPenalty the gap penalties used during alignment
128         * @param subMatrix the set of substitution scores used during alignment
129         */
130        protected AbstractProfileProfileAligner(Future<ProfilePair<S, C>> query, Profile<S, C> target,
131                        GapPenalty gapPenalty, SubstitutionMatrix<C> subMatrix) {
132                super(gapPenalty, subMatrix);
133                queryFuture = query;
134                this.target = target;
135                reset();
136        }
137
138        /**
139         * Sets the query {@link Profile}.
140         *
141         * @param query the first {@link Profile} of the pair to align
142         */
143        public void setQuery(Profile<S, C> query) {
144                this.query = query;
145                queryFuture = null;
146                reset();
147        }
148
149        /**
150         * Sets the target {@link Profile}.
151         *
152         * @param target the second {@link Profile} of the pair to align
153         */
154        public void setTarget(Profile<S, C> target) {
155                this.target = target;
156                targetFuture = null;
157                reset();
158        }
159
160        // method for ProfileProfileAligner
161
162        @Override
163        public ProfilePair<S, C> getPair() {
164                if (pair == null) {
165                        align();
166                }
167                return pair;
168        }
169
170        // methods for ProfileProfileScorer
171
172        @Override
173        public Profile<S, C> getQuery() {
174                return query;
175        }
176
177        @Override
178        public Profile<S, C> getTarget() {
179                return target;
180        }
181
182        // methods for AbstractMatrixAligner
183
184        @Override
185        protected CompoundSet<C> getCompoundSet() {
186                return (query == null) ? null : query.getCompoundSet();
187        }
188
189        @Override
190        protected List<C> getCompoundsOfQuery() {
191                // TODO replace with consensus sequence
192                return (query == null) ? new ArrayList<C>() : query.getAlignedSequence(1).getAsList();
193        }
194
195        @Override
196        protected List<C> getCompoundsOfTarget() {
197                // TODO replace with consensus sequence
198                return (target == null) ? new ArrayList<C>() : target.getAlignedSequence(1).getAsList();
199        }
200
201        @Override
202        protected int[] getScoreMatrixDimensions() {
203                return new int[] { query.getLength() + 1, target.getLength() + 1, (getGapPenalty().getType() == Type.LINEAR) ?
204                                1 : 3 };
205        }
206
207        @Override
208        protected int getSubstitutionScore(int queryColumn, int targetColumn) {
209                return getSubstitutionScore(qfrac[queryColumn - 1], tfrac[targetColumn - 1]);
210        }
211
212        @Override
213        protected boolean isReady() {
214                // TODO when added to ConcurrencyTools, log completions and exceptions instead of printing stack traces
215                try {
216                        if (query == null && queryFuture != null) {
217                                query = queryFuture.get();
218                        }
219                        if (target == null && targetFuture != null) {
220                                target = targetFuture.get();
221                        }
222                        reset();
223                } catch (InterruptedException e) {
224                        logger.error("Interrupted Exception: ", e);
225                } catch (ExecutionException e) {
226                        logger.error("Execution Exception: ", e);
227                }
228                return query != null && target != null && getGapPenalty() != null && getSubstitutionMatrix() != null &&
229                                query.getCompoundSet().equals(target.getCompoundSet());
230        }
231
232        @Override
233        protected void reset() {
234                super.reset();
235                pair = null;
236                if (query != null && target != null && getGapPenalty() != null && getSubstitutionMatrix() != null &&
237                                query.getCompoundSet().equals(target.getCompoundSet())) {
238                        int maxq = 0, maxt = 0;
239                        cslist = query.getCompoundSet().getAllCompounds();
240                        qfrac = new float[query.getLength()][];
241                        for (int i = 0; i < qfrac.length; i++) {
242                                qfrac[i] = query.getCompoundWeightsAt(i + 1, cslist);
243                                maxq += getSubstitutionScore(qfrac[i], qfrac[i]);
244                        }
245                        tfrac = new float[target.getLength()][];
246                        for (int i = 0; i < tfrac.length; i++) {
247                                tfrac[i] = target.getCompoundWeightsAt(i + 1, cslist);
248                                maxt += getSubstitutionScore(tfrac[i], tfrac[i]);
249                        }
250                        max = Math.max(maxq, maxt);
251                        score = min = isLocal() ? 0 : (int) (2 * getGapPenalty().getOpenPenalty() + (query.getLength() +
252                                        target.getLength()) * getGapPenalty().getExtensionPenalty());
253                }
254        }
255
256        // helper method that scores alignment of two column vectors
257        private int getSubstitutionScore(float[] qv, float[] tv) {
258                float score = 0.0f;
259                for (int q = 0; q < qv.length; q++) {
260                        if (qv[q] > 0.0f) {
261                                for (int t = 0; t < tv.length; t++) {
262                                        if (tv[t] > 0.0f) {
263                                                score += qv[q]*tv[t]*getSubstitutionMatrix().getValue(cslist.get(q), cslist.get(t));
264                                        }
265                                }
266                        }
267                }
268                return Math.round(score);
269        }
270
271}