Øvelse: BLAST
Kommenteret
version - marts 2008
Øvelse skrevet af: Rasmus
Wernersson
Introduktion
I denne øvelse skal vi bruge BLAST (Basic Local Alignment Search Tool) til at søge i
sekvensdatabaser såsom GenBank (DNA) og UniProt (protein). Et
meget væsentligt aspekt mht. brug af BLAST er at kunne vurderer signifikansen af et givet resultat.
BLAST er en program-pakke der kan hentes kvit og frit, og som
indeholder programmer til både at oprette nye
søge-database (fortmatdb
programmet) og til at søge i disse database (fx blastn til DNA-database og blastp til proteindatabaser).
Vi skal her bruge det webinterface til BLAST som findes hos NCBI: NCBI
BLAST - en af de store fordele er (ud over at vi ikke skal rode med en
UNIX prompt) at alle NCBI's databaser (+ UniProt/SwissProt !) er
tilgængelige til at søge i.
Links:
Del 1 - alignment af tilfældige sekvenser.
Som gennemgået i forelæsningen, vil er være risiko
for at få falske positive
resultater (altså et hit til en ikke-relateret sekvens) af rent
stokastiske årsager. I første del af øvelsen skal
vi kigge på hvad der sker hvis man søge med rent
tilfældige sekvenser som input.
I har alle fået udleveret en fire-sidet og en tyve-sidet
terning, som vi skal bruge til at generere DNA- og protein-sekvenser
(hvis du laver denne øvelser derhjemme så kan du istedet
bruge dette lille Java program "Remote
Roller") .
TRIN 1 - DNA sekvenser og BLASTN:
- Start med at slå 3 stk 25bp DNA sekvenser, ud fra
følgende tabel.

1 : A
2 : C
3 : G
4 : T
- Vi skal nu have fat i BLASTN algoritmen hos NCBI.
- Vælg "Nucleotide BLAST"
fra BLAST hovedsiden.
- Under "Program Selection"
skal I vælge "Somewhat similar
sequences (blastn)"
- Vælg "Nucletide
Collection NR/NT" som database. NR står for "Non
Redundancy" og er en samling af alle
ikke-redundante sekvenser fra GenBank og genom-sekventerings
projekterne.
- VIGTIGT: Til denne kørsel med de små kunstige DNA
sekvenser, skal vi slå noget af NCBI automatik mht. valg af
options fra. Tryk på "Algorithm
parameters"for at få adgang til de bagvedliggende options
(se nedenstående screenshot). Fravælg "Automaticlly adjust parameters for short
input sequences". Bemærk:
her kan man også indstille E-value cutoff ("Expect threshold") og BLASTN "word
size".
- Indsæt dine tre
sekvenser (hint: FASTA format) - og start BLAST søgningen.

Svar:
Husk
at det er nemt at skrive sekvenser op i FASTA format i hånden -
hver sekvens skal blot have en header der fortæller om dens navn:
>seq_01
CGCCCGACCGTGTAGGAGGTCCGGT
>seq_02
TTCCAAGTAGGTTATGTTAAAGCCG
>seq_03
CGGTTAATTCGCTGTCTATATAAGG
- Lad os kigge lidt på resultaterne - bemærk at der
både er en tabel med en oversigt i toppen, samt at selve de
individuelle alignments vises længere nede på siden:
- Hvor stor (i basepar) er den database vi har BLAST'et imod?
Svar: 23
milliarder bogstaver (dvs. baser i tilfældet med DNA).
- Er der nogen sekvenser de ligner vores tilfældigt
genererede?
Svar:
Der er flere perfekte hits, fx:
GENE ID: 5438907 BC1G_02961 | hypothetical
protein
[Botryotinia fuckeliana B05.10]
Score = 33.7 bits (36), Expect = 4.8
Identities = 18/18 (100%), Gaps = 0/18 (0%)
Strand=Plus/Minus
Query 1 CGCCCGACCGTGTAGGAG 18
||||||||||||||||||
Sbjct 574 CGCCCGACCGTGTAGGAG 557
- Hvor lange er de hits der findes? Svar: Typisk en
17-19 baser.
- Hvor mange % identity? Svar: 95% - 100%
- I hvilket område ligger bit-scores ("max score")? Svar: typisk
på 33 - 35 bits.
- I hvilket område ligger E-value?
Svar:
1.4 - 4.8 (for denne søgning, ligger typisk fra 0.1 - 10.0).
BEMÆRK: at default cut-off for e-value er 10.0 - vi ser derfor
ikke hits der har en dårligere (større) e-value end 10.0.
- Giver disse resultater nogen biologisk
mening?
Svar:
NEJ! - hits'ne er rigtige nok, det er faktisk sekvens der findes
i databasen. Men vi ved jo at vores sekvenser er helt tilfældige,
og derfor ikke har noget slægtskab med dem i databasen. Vi har
kun fundet vores hits, fordi databasen er så enorm stor, at vi af
rent stokastiske årsager rammer nogle sekvenser der ligner.
E-værdierne fortæller netop dette: Som jeg gennemgik i
forelæsningenen, vil alignment
score ("bit score"
i BLAST) følge en ekstremværdi-fordeling for de sekvenser
der ikke har noget med vores input-sekvenser at gøre (fordelt
omkring en lav bit-score), og de hits der faktisk er relateret
følger en normalfordeling (omkring en væsentlig
højere bit-score).
Ud fra ekstremværdi-fordelingen (som beregnes når BLAST
databasen oprettes), er det mulig at estimere hvor mange
tilfældige, ikke-relaterede hits, man må forvente ved en
given bit-score. Dette udtrykkes som en e-værdi (EXPECT-value)
som fortæller, hvor mange ikke-relaterede hits man kan forvente at finde ved en given
bit-score. Hvis vi fx kigger på vores ovenstående hit med
en bit-score på 33.7, giver dette en e-value på 4.8 (i
denne database, af 23*109 baser). Man må altså
kunne forvente statistisk set af finde 4.8 hits med en score på
33.7 bits i denne database - af
tilfældige årsager!
- Lad os prøve at vælgde en anden database. Vælg
"Human genomic plus transcripts"
som database og kør søgningen en gang til.
- En kort kommentar om word-size
og DNA sekvenser: Som default er word-size 11 for DNA sekvenser - dette
betyder at algorithmen først kigger efter et perfect match på 11 baser, og
hvis det findes, så kigges der på hvad der ligger af
sekvens omkring dette stykke - man kan ændre den til 7
(langsommere, men mere grundigt) eller 15 (hurtigere - finder ikke helt
så meget). Det gør ikke den store forskel i dette eksempel
med tilfældige sekvenser, da de samle hits altid er korte og
næsten altid uden gaps - men det kan betyde en forskel, hvis man
aligner rigtigt DNA, hvor resultatet er et længe alignment med
både mismatches og gaps.
Kommentar:
Som gennemgået i forelæsningen - og i PodCast'en.
TRIN 2 - Protein sekvenser og BLASTP:
- Lad os nu istedet kaste os over nogle protein-sekvenser. Start
med at slå 3 peptider med en længe på 25 aa, ud fra
nedenstående tabel.

1
: A
2 : R
3 : N
4 : D
5 : C
6 : Q
7 : E
8 : G
9 : H
10 : I
11 : L
12 : K
13 : M
14 : F
15 : P
16 : S
17 : T
18 : W
19 : Y
20 : V
- Åbn "Protein BLAST"
siden hos NCBI og vælg BLASTP som metode. Sæt sekvenserne
ind og vælg igen "NR" databasen (denne gang består den af
oversætte CDS'er + SwissProt mm).
Under "Algorithm Parameter"vælges
BLOSUM62 som matrice, og "Expect treshold" sættes til 1000 (default:
10) VIGTIGT. Husk også at slå "Short sequence"
autoparameterne fra (ligesom med DNA søgningen) - ellers
ignoreres vores omhyggeligt indstillede parametere.
Submit derefter søgningen.
HUSK E-VALUE
CUT-OFF PÅ 1000!
Mine
sekvenser:
>seq_01
aymgsplsflshdhncirdakyvvs
>seq_02
danshrsqwnrlgmanwfqhwsewv
>seq_03
vsltdgpdkhesqhirqdkyfddia
- Hvor stor er databasen? Svar: 2.1 milliader
bogstaver (aminosyrere).
- Hvor lange er hits'ne typisk og har de gaps? Svar: Typisk 17-21.
Ingen gaps, men med en del mismatches.
- Hvor mange % identity? Svar: ca. 40-60
- I hvilket område ligger E-value? Svar: 229 - 986 (!).
- Prøve at kigge grundigt på et par af selve
alignments'ne - ("+" betyder similar) - er er nogle af dem der som
sådan ser rimeligt plausible ud (hvis vi for et øjeblik
ser bort fra længden)?
Svar: Ja
se fx. følgende alignment - det ser da overbevisnede ud
(men det er alt, alt for kort til at være signifikant):
Score
= 27.7 bits (60), Expect = 174, Method:
Composition-based stats.
Identities = 10/17 (58%),
Positives = 12/17 (70%), Gaps = 0/17 (0%)
Query 1
AYMGSPLSFLSHDHNCI 17
AY GSP + HDHNC+
Sbjct 49
AYPGSPHGYDIHDHNCL 65
- Hvis vi havde brugt default cut-off'et for E-value på 10,
ville der så være fundet nogen sekvenser?
Svar:
NEJ. Hvis vi også havde kørt søgningerne mod
DNA databasen med en e-value cut-off på 1000, ville vi har set
sidevis af hits for hver sekvens.
- Fanger vi mest junk ved at bruge DNA eller peptid i vores
søgninger?
Svar:
Risikoen for at samle et "falsk" hit op, er større med DNA
- bemærk at vi var nødt til at "snyde" med cut-off for i
det hele taget at se nogle hits mod proteinsekvenserne.
- Kommentar ang. word-size
og proteinsekvenser. blastp
fungerer noget anderledes end blastn
mht. at finde et sted et sekvenserne med høj-similariet: istedet
for at kræve et længede stykke der er 100% det samme (det
vil ikke give mening for protein-sekvens), prøver algoritmen at
finde to sekvens-stumper (af længde "word-size", default =3) med
høj similaritet inden for en kort afstand (default: 40 aa).
Så snart et sådan sekvens par er fundet, fylder BLAST ind
med almindeligt lokalt alignment.
Kommentar:
Som gennemgået i forelæsningen - og i PodCast'en.
Del 2 - overførelse af viden om funktion
En af de allermest brugte anvendelser af BLAST er til, at
overføre viden om funktion fra en sekvens til en anden. En
vigtig forudsætning for dette er erkendelsen af alle levende
organismer er evolutionært beslægtede. Hvis gen "X" i en
organisme fx. vides at kode for et bestemt enzym (e.g. en protease) kan
man med rimelig sikkerhed antage at gen "Y" fra en anden organisme har
sammen funktion, hvis de to gener er tilstrækkeligt ens (hvilket
betyder at de ikke er gået alt for meget deres egne veje, siden
sidste fælles stamform for genet).
Vi undersøgte præcis dette fænomen i "Parvis
alignment" øvelsen (med en serin-protease fra to Bacills arter). Denne gang har vi
ikke to sekvenser vi iforvejen har en formodning om er beslægtede
- vi skal derimod bruge BLAST til at finde de sekvenser der ligner i de
samlede database med millioner af sekvenser.
Lad os starte med en situation der giver nogle ganske fornuftige
hits i databasen. Nedenstående sekvens er et komplet transcript
(mRNA) fra en bakterie. Lad os finde ud af hvad det er.
>Unknown_transcript01
CCACTTGAAACCGTTTTAATCAAAAACGAAGTTGAGAAGATTCAGTCAACTTAACGTTAATATTTGTTTC
CCAATAGGCAAATCTTTCTAACTTTGATACGTTTAAACTACCAGCTTGGACAAGTTGGTATAAAAATGAG
GAGGGAACCGAATGAAGAAACCGTTGGGGAAAATTGTCGCAAGCACCGCACTACTCATTTCTGTTGCTTT
TAGTTCATCGATCGCATCGGCTGCTGAAGAAGCAAAAGAAAAATATTTAATTGGCTTTAATGAGCAGGAA
GCTGTTAGTGAGTTTGTAGAACAAGTAGAGGCAAATGACGAGGTCGCCATTCTCTCTGAGGAAGAGGAAG
TCGAAATTGAATTGCTTCATGAATTTGAAACGATTCCTGTTTTATCCGTTGAGTTAAGCCCAGAAGATGT
GGACGCGCTTGAACTCGATCCAGCGATTTCTTATATTGAAGAGGATGCAGAAGTAACGACAATGGCGCAA
TCAGTGCCATGGGGAATTAGCCGTGTGCAAGCCCCAGCTGCCCATAACCGTGGATTGACAGGTTCTGGTG
TAAAAGTTGCTGTCCTCGATACAGGTATTTCCACTCATCCAGACTTAAATATTCGTGGTGGCGCTAGCTT
TGTACCAGGGGAACCATCCACTCAAGATGGGAATGGGCATGGCACGCATGTGGCCGGGACGATTGCTGCT
TTAAACAATTCGATTGGCGTTCTTGGCGTAGCGCCGAGCGCGGAACTATACGCTGTTAAAGTATTAGGGG
CGAGCGGTTCAGGTTCGGTCAGCTCGATTGCCCAAGGATTGGAATGGGCAGGGAACAATGGCATGCACGT
TGCTAATTTGAGTTTAGGAAGCCCTTCGCCAAGTGCCACACTTGAGCAAGCTGTTAATAGCGCGACTTCT
AGAGGGGTTCTTGTTGTAGCGGCATCTGGGAATTCAGGTGCAGGCTCAATCAGCTATCCGGCCCGTTATG
CGAACGCAATGGCAGTCGGAGCGACTGACCAAAACAACAACCGCGCCAGCTTTTCACAGTATGGCGCAGG
GCTTGACATTGTCGCACCAGGTGTAAACGTGCAGAGCACATACCCAGGTTCAACGTATGCCAGCTTAAAC
GGTACATCGATGGCTACTCCTCATGTTGCAGGTGCAGCAGCCCTTGTTAAACAAAAGAACCCATCTTGGT
CCAATGTACAAATCCGCAATCATCTAAAGAATACGGCAACGAGCTTAGGAAGCACGAACTTGTATGGAAG
CGGACTTGTCAATGCAGAAGCGGCAACACGCTAATCAATAATAATAGGAGCTGTCCCAAAAGGTCATAGA
TAAATGACCTTTTGGGGTGGCTTTTTTACATTTGGATAAAAAAGCACAAAAAAATCGCCTCATCGTTTAA
AATGAAGGTACC
- Søg efter sekvensen i NR/NT
databasen (BLASTN) med default
indstillinger.
- Lad os prøve det samme med sekvensen på
protein-niveau.
- Oversæt den længste ORF med VirtualRibosome
(tip: husk at vælge alle
læserammer), og gem/kopier
proteinsekvesen i FASTA format.
Svar:
HUSK at bruger OFR finderen!
Eftersom vi ved at det er et fuld-længde transcript, kan vi godt
regne med at både START og STOP codon er med. Vi kan derfor
vælge "START-codon = strict".
Bemærk:
man kan oversætte med både den standard-genetiske kode
(tabel 1) eller alternativt med tabel 11 ("The Bacterial and Plant
Plastid Code"). Den eneste forskel er at tabel 11 tillader nogle ekstra
yderst sjældne START codons, bla. fra bacterielle phager.
- Brug BLASTP mod NR.
- Findes der nogen konserverede protein-domæner? Dette kan
være en vigtig brik i puslespillet med at finde ud af hvad en
given proteinsekvens har af funktion.
Svar:
Ja, der findes et "Peptidase S8" domæne. Dette kan se på
både "pause-siden" mens der søges, og ved at klikke
på "Show conserved domains"
øverst på resultat siden.
- Er der nogen signifikante hits?
Svar: Ja, ganske mange.
- Hvilken type enzymer er de bedste hits?
Svar:
Proteaser og Peptidaser (hvilket dækker over det samme).
Bemærk at man kan klikke på hits'ne og gå direkte til
sekvensen, hvor man kan læse yderligere. Vi kan slå
følgende fast: 1) Der er en serin protease, 2) af S8 familien,
3) den er basisk, 4) den kommer fra Bacillus.
- Hvordan ser det ud med fordelingen af hits (høj ->
lav kvalitet)?
Svar:
Fra fuld-længde hits med en e-value på mindre end 1e-169 og
gradvis
dernedaf (per default ser vi kun de 100 bedste hits).
- Er protein-alignment'et bedre til at fange breden af
similaritet? Svar:
Ja.
- OK - det var jo lidt snyd at bruge en sekvens, der allerede var
tilstede i databasen, så lad os istedet kaste os over
følgende sekvens.
Det er her tale om et DNA fragment sekventeret direkte fra DNA
ekstraheret fra en jordprøve. Fragmentet der går under det
poetiske navn "KLON12" er amplificeret med PCR primere der er
særligt designet til at amplificere kernen af gener der koder for en
bestem klasse af enzymer. (I parantes bemærket er det er
gen-fragment jeg selv har sekventeret - mit speciale handlede om
ikke-kultiverbare mikroorganismer fra jord).
LOCUS
KLON12.DNA 609 BP
DS-DNA
UPDATED 06/14/98
DEFINITION UWGCG file capture
ACCESSION -
KEYWORDS -
SOURCE -
COMMENT
Non-sequence data from original file:
BASE
COUNT 174 A 116
C 162 G 157
T 0 OTHER
ORIGIN ?
klon12.dna Length:
609 Jun 13, 1998 - 03:39 PM Check: 6014 ..
1 AACGGGCACG GGACGCATGT AGCTGGAACA GTGGCAGCCG TAAATAATAA TGGTATCGGA
61 GTTGCCGGGG TTGCAGGAGG AAACGGCTCT ACCAATAGTG GAGCAAGGTT AATGTCCACA
121
CAAATTTTTA ATAGTGATGG GGATTATACA AATAGCGAAA CTCTTGTGTA CAGAGCCATT
181
GTTTATGGTG CAGATAACGG AGCTGTGATC TCGCAAAATA GCTGGGGTAG TCAGTCTCTG
241
ACTATTAAGG AGTTGCAGAA AGCTGCGATC GACTATTTCA TTGATTATGC AGGAATGGAC
301
GAAACAGGAG AAATACAGAC AGGCCCTATG AGGGGAGGTA TATTTATAGC TGCCGCCGGA
361
AACGATAACG TTTCCACTCC AAATATGCCT TCAGCTTATG AACGGGTTTT AGCTGTGGCC
421
TCAATGGGAC CAGATTTTAC TAAGGCAAGC TATAGCACTT TTGGAACATG GACTGATATT
481
ACTGCTCCTG GCGGAGATAT TGACAAATTT GATTTGTCAG AATACGGAGT TCTCAGCACT
541
TATGCCGATA ATTATTATGC TTATGGAGAG GGAACATCCA TGGCTTGTCC ACATGTCGCC
601
GGCGCCGCC
//
- Opgaven er nu at finde
ud af hvad denne stump sekvens koder for,
ved at bruge de metoder vi netop har gennemgået.
Bemærk: Sekvensen er
(mere eller mindre) i GenBank format og NCBI BLAST forventer enten
FASTA eller rå sekvens. De fleste webservere (ikke alle!), der
også kan acceptere en sekvens i rå format, sørger
automatisk for at smide alle mellemrum og tal væk. Dette
gælder også for både NCBI BLAST og de fleste af de
servere jeg selv har skrevet, fx VirtualRibosome. Man kan altså
kopiere de linier ind hvor selve sekvensen står. NB NB NB: Hvis man skal submitte
flere sekvenser på een gang, virker dette trick ikke - brug
korrekt formaterede FASTA sekvenser istedet.
Tænk over følgende inden I går i gang med
søgningen:
- Er der tale om proteinkodende sekvens? Svar: Ja.
- Hvis ja, kan man så regne med at START og STOP codon er
med?
Svar:
Nej - det er et fragment af en gen - det ved vi kun fordi jeg
allerede have giver jer en hint om det. Hvis man ikke har den oplysning
kan man prøve sig frem.
Lad os som udgangspunkt sige at vi ikke tror på noget med en
E-værdi på større end ("dårligere end") 10-10.
Dette kunne være et typisk eksempel på en eksamensopgave.
Det er vigtigt at I i jeres svar argumenterer (bare i stikords form,
for jeres egen skyld, siden I ikke skal aflevere denne opgave) for
jeres strategi for at nå frem til et resultat:
- Hvilke værktøjer / databaser er relevante at bruge?
Svar: Vi
vil i sidste ende gerne bruge BLAST til at søge i de store
samlede databaser.
Lad os derfor prøve følgende:
1) BLASTN
2) Oversætte til protein (fx. VirtualRibosome).
3) BLASTP
Både med BLASTN og BLASTP ønsker vi at bruge NR databasen
- det er den der giver den største bredde. Der er godt at
started med en stor, bred database - så kan man evt. gå i
dybten i en af de mere specifikke, når man har skudt sig ind
på hvad for en type sekvens man har.
- Hvad er resultaterne af de forskellige typer søgning -
hvad virker og hvad virker ikke (og hvorfor)?
Svar:
1) BLASTN - når vi prøver at BLAST'e mod NR får vi
slet ikke nogen signifikante hits:
gb|CP000159.1|
Salinibacter ruber DSM 13855, complete
genome 51.8 0.003
gb|CP000813.1| Bacillus pumilus SAFR-032, complete
genome
50.0 0.012
gb|CP000485.1| Bacillus thuringiensis str. Al Hakam,
complete... 48.2 0.042
gb|CP000001.1| Bacillus cereus E33L, complete
genome
48.2 0.042
dbj|AB085752.1| Bacillus sp. KSM-LD1 gene for protease,
compl... 48.2 0.042
ref|XM_001193501.1| PREDICTED: Strongylocentrotus purpuratus
... 44.6 0.51
ref|XM_001196149.1| PREDICTED: Strongylocentrotus purpuratus
... 44.6 0.51
gb|BT014524.1| Lycopersicon esculentum clone 133938F, mRNA
se... 44.6 0.51
emb|AM486171.2| Vitis vinifera contig VV78X025776.9, whole
ge... 42.8 1.8
emb|AM478115.2| Vitis vinifera contig VV78X115679.5, whole
ge... 42.8 1.8
emb|CR392038.9| Zebrafish DNA sequence from clone DKEY-172J4
... 42.8 1.8
ref|XM_001191953.1| PREDICTED: Strongylocentrotus purpuratus
... 42.8 1.8
dbj|AB214317.1| Bacillus licheniformis prz gene for
keratinaz... 42.8 1.8
gb|AY136509.1| Glycine max subtilisin-like protease C1 gene,
... 42.8 1.8
gb|AC133779.21| Medicago truncatula clone mth2-27n4,
complete... 42.8 1.8
gb|AF036960.3| Glycine max subtilisin-like protease C1 mRNA,
... 42.8 1.8
ref|XM_779662.1| PREDICTED: Strongylocentrotus purpuratus
sim... 42.8 1.8
emb|CU329670.1| Schizosaccharomyces pombe chromosome
I
41.0 6.2
dbj|AB305158.1| Alkalimonas collagenimarina aprI gene for
sub... 41.0 6.2
emb|CT030005.14| Zebrafish DNA sequence from clone
DKEY-206D1... 41.0 6.2
ref|XM_001583403.1| Trichomonas vaginalis G3 Dynein heavy
cha... 41.0 6.2
ref|XM_001180101.1| PREDICTED: Strongylocentrotus purpuratus
... 41.0 6.2
ref|XM_785195.2| PREDICTED: Strongylocentrotus purpuratus
sim... 41.0 6.2
ref|XR_023110.1| PREDICTED: Pan troglodytes kinase
suppressor... 41.0 6.2
ref|XM_001105629.1| PREDICTED: Macaca mulatta kinase
suppress... 41.0 6.2
ref|NM_115638.1| Arabidopsis thaliana unknown protein
(AT3G57... 41.0 6.2
dbj|AK124344.1| Homo sapiens cDNA FLJ42353 fis, clone
UTERU20... 41.0 6.2
ref|XM_388640.1| Gibberella zeae PH-1 strain PH-1; NRRL
31084... 41.0 6.2
gb|AC108679.6| Homo sapiens 3 BAC RP11-485G4 (Roswell Park
Ca... 41.0 6.2
gb|AC148994.3| Medicago truncatula chromosome 7 clone
mth2-3c... 41.0 6.2
gb|AY126451.1| Dactylaria parvispora serine protease mRNA,
pa... 41.0 6.2
(TIP:
Man kan få dette tekst-baserede output ved at klikke på "Reformat these result" i toppen af
siden, og så vælge "Plain text").
Der er simpelt hen ikke noget i databasen, der er ligner vores input
tilstrækkeligt - derfor finder vi kun de tilfældige
resultater. En søgning på DNA niveau er kun egnet til at
finde de sekvenser, der ligner virkeligt meget.
2) Virtual ribosome - resultat fra ORF finderen (NB: tabel 1, "Standard
Genetic Code" kan også sagtens bruges).
VIRTUAL RIBOSOME
----------------
Translation table: Bacterial and Plant
Plastid
>Seq1_rframe1_ORF
Reading frame: 1
N G
H G T H V A G T V
A A V N N N G I G
V A G V A G G N G
S
5'
AACGGGCACGGGACGCATGTAGCTGGAACAGTGGCAGCCGTAAATAATAATGGTATCGGAGTTGCCGGGGTTGCAGGAGGAAACGGCTCT
90
..............................))).....................))).................................
T N
S G A R L M S T Q
I F N S D G D Y T
N S E T L V Y R A
I
5'
ACCAATAGTGGAGCAAGGTTAATGTCCACACAAATTTTTAATAGTGATGGGGATTATACAAATAGCGAAACTCTTGTGTACAGAGCCATT
180
.....................>>>.........))).......................................))).........)))
V Y
G A D N G A V I S
Q N S W G S Q S L
T I K E L Q K A A
I
5'
GTTTATGGTGCAGATAACGGAGCTGTGATCTCGCAAAATAGCTGGGGTAGTCAGTCTCTGACTATTAAGGAGTTGCAGAAAGCTGCGATC
270
........................))))))...........................)))...)))......)))............)))
D Y
F I D Y A G M D E
T G E I Q T G P M
R G G I F I A A A
G
5'
GACTATTTCATTGATTATGCAGGAATGGACGAAACAGGAGAAATACAGACAGGCCCTATGAGGGGAGGTATATTTATAGCTGCCGCCGGA
360
.........)))............>>>...............)))............>>>.........)))...)))............
N D
N V S T P N M P S
A Y E R V L A V A
S M G P D F T K A
S
5'
AACGATAACGTTTCCACTCCAAATATGCCTTCAGCTTATGAACGGGTTTTAGCTGTGGCCTCAATGGGACCAGATTTTACTAAGGCAAGC
450
........................>>>...........................)))......>>>........................
Y S
T F G T W T D I T
A P G G D I D K F
D L S E Y G V L S
T
5'
TATAGCACTTTTGGAACATGGACTGATATTACTGCTCCTGGCGGAGATATTGACAAATTTGATTTGTCAGAATACGGAGTTCTCAGCACT
540
...........................)))..................)))............)))........................
Y A
D N Y Y A Y G E G
T S M A C P H V A
G A A
5'
TATGCCGATAATTATTATGCTTATGGAGAGGGAACATCCATGGCTTGTCCACATGTCGCCGGCGCCGCC
609
.......................................>>>...........................
(TIP:
Husk at du kan få sekvensen i FASTA format ved at klikke på
"FASTA" linket på resultat-siden).
3) BLASTP
Det
første vi ser, er at at der finds et "Peptidase_S8"
domæne. Dette er en meget kraftig indikator på funktionen.
Vi får en række signifinkante hits. Hvis vi kigger på
de bedste, ser at de alle er serin-proteaser (af S8 familien).
ref|NP_809125.1|
protease [Bacteroides thetaiotaomicron VPI-5...
154 4e-36
ref|ZP_02425799.1| hypothetical protein ALIPUT_01953
[Alistip... 147 5e-34
ref|ZP_02437453.1| hypothetical protein BACSTE_03728
[Bactero... 142 2e-32
ref|ZP_02424566.1| hypothetical protein ALIPUT_00684
[Alistip... 137 3e-31
ref|YP_446403.1| peptidase families S8 and S53 domain
protein... 129 1e-28
ref|ZP_02063425.1| hypothetical protein BACOVA_00373
[Bactero... 126 7e-28
gb|AAB49694.1| protease [Bacillus sp.
PD498]
117 4e-25
ref|YP_475789.1| peptidase, S8A (subtilisin) family
[Synechoc... 111 3e-23
ref|YP_477062.1| peptidase, S8A (subtilisin) family
[Synechoc... 110 6e-23
ref|NP_693519.1| thermophilic serine proteinase
[Oceanobacill... 107 5e-22
gb|AAB36499.1| thermostable alkaline protease
[Thermoactinomy... 103 5e-21
sp|Q45670|THES_BACSJ Thermophilic serine proteinase
precursor... 103 7e-21
gb|AAK29176.1| thermophilic alkaline protease [Geobacillus
st... 103 8e-21
ref|ZP_01860436.1| thermophilic serine proteinase [Bacillus
s... 103 9e-21
pdb|1DBI|A Chain A, Crystal Structure Of A Thermostable
Serin... 102 2e-20
ref|ZP_01170337.1| subtilisin-type proteinase [Bacillus sp.
N... 101 3e-20
ref|ZP_02328516.1| subtilisin-type proteinase [Paenibacillus
... 100 5e-20
ref|NP_487583.1| protease [Nostoc sp. PCC 7120]
>dbj|BAB75242... 100 5e-20
ref|ZP_00106975.2| COG1404: Subtilisin-like serine proteases
... 100 5e-20
ref|YP_324024.1| Peptidase S8 and S53, subtilisin, kexin,
sed... 100 6e-20
ref|ZP_01630730.1| protease [Nodularia spumigena CCY9414]
>gb... 100 8e-20
ref|YP_175254.1| serine protease [Bacillus clausii KSM-K16]
>... 99.8 9e-20
ref|NP_242357.1| protease [Bacillus halodurans C-125]
>dbj|BA... 99.8 9e-20
(BEMÆRK:
"Hypotetical" proteiner er ikke nær så stærk evidens
som et velbeskrevet protein).
Lad os kigge lidt nærme på det bedste hit (NP_809125):
>ref|NP_809125.1|
protease [Bacteroides thetaiotaomicron VPI-5482]
gb|AAO75319.1| protease [Bacteroides thetaiotaomicron VPI-5482]
Length=695
Score = 154 bits (389), Expect = 4e-36, Method:
Compositional matrix adjust.
Identities = 98/238 (41%), Positives = 134/238 (56%), Gaps =
45/238 (18%)
Query 2
GHGTHVAGTVAAVNNNGIGVAGVAGGNGSTNSGARLMSTQIFNSDGDYT-NSETLVYRAI 60
GHGTHVAGT+AA+NNNG GV G+AGG+GS NSG ++MS Q+F + T
++E RAI
Sbjct 295
GHGTHVAGTIAAMNNNGEGVCGIAGGDGSKNSGVKIMSCQVFAGEAGVTLDAEA---RAI 351
Query 61
VYGADNGAVISQNSWGSQS-------------------LTIKELQKAAIDYFIDYAGMDE 101
Y ADNGAVI Q SWG
S
+ L+K A+DYFI+ AG
Sbjct 352
KYAADNGAVILQCSWGYNSSLANLIEGYTPGPGSEEEWEKLYPLEKDALDYFINNAGS-- 409
Query 102
TGEIQTGPMRGGIFIAAAGNDNVSTPNMPSAYERVLAVASMGPDFTKASYSTFGTWTDIT 161
G + GG+ I A+GN+ P+AY + ++V+++
DFT ASYS +G I+
Sbjct 410
----PNGVIDGGLAIFASGNEYAGMAAFPAAYSKCISVSAVAADFTPASYSNYGKEVTIS 465
Query 162
APGGDIDKFD-------------LSEYGVLSTYADN---YYAYGEGTSMACPHVAGAA 203
APGGD +
++
+ +LST+ N Y + +GTSMACPHV+G A
Sbjct 466
APGGDTEYYNPVGQDDPEGWEEGIHSGSILSTWIQNGNATYGFMDGTSMACPHVSGVA 523
Bemærk at selv om det ikke er et perfekt hit (vores input-sekvens
findes jo netop ikke i databasen), så ser det ganske fornuftigt
ud: alignment'et er langt med en Identity på 41% og en Similarity
på 56%. Sammenholdt med det faktum at alle de bedste hits er
serin-proteaser, har vi en meget stærk indikation på
hvilken type af enzym vi har med at gøre.
- Endelig konklusion: Hvad er KLON12?
Svar: En
serin-protease, fra S8 familien.
Del 3 - BLAST og genomer
Nu har vi indtil videre brugt BLAST til at søge i de store brede
databaser (GenBank + andre DNA-databaser, UniProt + andre
proteindatabaser). Lad os runde øvelsen af med at kigge på
hvordan man bruger BLAST på nogle mere smalle databaser, nemlig
til at søge i specifikke genomdatabaser.
Dette vil typisk være nyttigt hvis man har en sekvens med kendt
funktion fra en bestemt organisme (fx. et gen der styrer cell-deling
hos gær, Saccharomyces
cerevisiae), og gerne vil finde de tilsvarende gener hos en
anden organisme (fx. menneske - gener der styrer cell-deling er ofte
indvolverede i kræft).
Som I sikkert har lagt mærke til har det hos NCBI været
muligt at søge i både musen og menneskets genom med BLAST.
Der har også været mulighed for kun at ville se resultater
fra en bestemt organisme.
Der findes en række andre servere der specifikt handler om en
bestemt organismes genom. For eksempel har gær-genomet "hjemme" i
Saccharomyces Genome Database
(SGD - www.yeastgenome.org)
- bemærk iøvrigt at også denne side tilbyder BLAST
som en måde at søge i databasen på. BLAST er ganske
enkelt standard-værktøjet til at finde sekvenser ud fra
andre sekvenser.
Vi skal her nøjes med at kigge på de genom-ressourcer der
ligger hos NCBI (med en enkelt afstikker til SGD):
http://www.ncbi.nlm.nih.gov/Genomes/
- Lad os undersøge hvordan histonerne (ansvarlige for at
pakke DNA) er beslægtede mellem gær og menneske
(evolutionær afstand: 1-1,5 milliard år).
Start med at slå HTA2
genet op i SGD (hint: der er en Quick Search box øverst på
siden). Bemærk at der er en kort forklaring om hvad genet og dets
protein-produkt gør. Hvad står der her om dets forhold til
genet "HTA1"?
Svar: HTA2
og HTA1 er næsten ens ("nearly identical").
Slå nu protein-sekvensen op (se under "Protein Info &
Structure" ude til højre). Hold
vinduet åbent eller gem
sekvensen, vi skal bruge den om lidt.
- Gå nu til NCBI's genom side. Læg mærke til at
der ude til højre findes en liste af organismer (under "Organism-specific") hvor de
informationer der findes er vist med bogstavkode ("B" = BLAST). Hvis man her
vælger BLAST for en specifik organisme, får man (ikke
overraskende) mulighed for at søge i dette specifikke genom (og
de oversatte proteinsekvenser) med BLAST.
Inden vi søger i det humane genom, så lad os
undersøge hvad vi kan finde hvis vi søget i gærs
eget genom. Vælg derfor "B"
under "Yeast".
- Vælg "RefSeq protein" som database.
- Vælg "BLASTP" som program.
- Søg efter HTA2 protein-sekvensen.
- Hvor mange pålidelige hits findes der? Svar: 3
- Giver det mening ud fra det der stod under beskrivelsen af HTA2?
Svar: Ja
- HTA1 og HTA2 er meget ens (kun 2 aminosyrer er forskellige).
Som bonus finder vi også "H2AZ" - en anden variant af HTA2 der
noget mere divergeret i sekvens.
- Næste trin er at søge i det (oversatte) humane
genom.
Gå ind på BLAST siden for det humane genom.
- Vælg "RefSeq protein" som database.
- Bemærk at der er flere muligheder for valg af database
end hos gær - det skyldes simpelthen at man har dårligere
styr på det humane genom, og at man som bruger kan have brug at
at vælge aktivt mellem forskellige fortolkninger.
- Vælg "BLASTP" som program og søg efter homologer
til HTA2 hos mennesket.
HITS:
ref|NP_003503.1|
H2A histone family, member L [Homo
sapiens] 191
1e-49
ref|NP_003507.1|
H2A histone family, member O [Homo sapiens] ...
190 2e-49
ref|NP_254280.1|
histone H2a [Homo
sapiens]
190 3e-49
ref|NP_003500.1|
H2A histone family, member C [Homo sapiens] ...
190 3e-49
ref|NP_542163.1|
H2A histone family member [Homo
sapiens]
190 3e-49
ref|NP_734466.1|
histone H2A [Homo
sapiens]
189 4e-49
ref|NP_778235.1|
histone H2A [Homo
sapiens]
189 5e-49
ref|NP_002096.1|
H2A histone family, member X [Homo
sapiens] 189
6e-49
ref|NP_003508.1|
H2A histone family, member Q [Homo
sapiens] 189
6e-49
ref|NP_066409.1|
histone 1, H2ad [Homo
sapiens]
189 6e-49
ref|NP_066390.1|
H2A histone family, member A [Homo sapiens] ...
189 6e-49
ref|NP_066544.1|
H2A histone family, member E [Homo
sapiens] 188
8e-49
ref|NP_808760.1|
H2A histone family, member J isoform 2 [Homo sa
188 1e-48
ref|NP_060737.1|
H2A histone family, member J isoform 1 [Homo sa
184 2e-47
ref|NP_061119.1|
core histone macroH2A2.2 [Homo
sapiens]
150 2e-37
ref|NP_613258.2|
H2A histone family, member Y isoform 3 [Homo sa
149 7e-37
ref|NP_613075.1|
H2A histone family, member Y isoform 1 [Homo sa
149 7e-37
ref|NP_004884.1|
H2A histone family, member Y isoform 2 [Homo...
149 7e-37
ref|NP_036544.1|
H2A histone family, member V isoform 1 [Homo sa
130 3e-31
ref|NP_002097.1|
H2A histone family, member Z [Homo
sapiens] 129
6e-31
ref|NP_619541.1|
H2A histone family, member V isoform 2 [Homo sa
111 1e-25
ref|NP_958844.1|
H2A histone family, member V isoform 3 [Homo sa
108 1e-24
ref|XP_946222.1|
PREDICTED: similar to Histone H2a [Homo sapiens
104 2e-23
ref|NP_542451.1|
H2A histone family, member B3 [Homo sapiens]...
87.4 3e-18
ref|NP_001017990.1|
H2A histone family, member B1 [Homo sapiens]
84.7 2e-17
ref|XP_934613.2|
PREDICTED: similar to testis-specific histon...
77.4 3e-15
ref|NP_958925.1|
H2A histone family, member V isoform 5 [Homo sa
77.4 3e-15
ref|XP_944049.2|
PREDICTED: similar to H2A histone family, me...
65.5 1e-11
ref|XP_370909.1|
PREDICTED: similar to H2A histone family, me...
65.5 1e-11
ref|NP_958924.1|
H2A histone family, member V isoform 4 [Homo sa
63.9 3e-11
ref|NP_005624.2|
son of sevenless homolog 1 [Homo
sapiens]
47.8 2e-06
ref|NP_001018082.1|
BTB (POZ) domain containing 11 isoform 3 [Ho
42.7 7e-05
ref|NP_665803.1|
ankyrin repeat and BTB (POZ) domain containing
37.0 0.004
- Hvor mange pålidelige hit findes der?
Svar:
Lad os ligesom før kræve at hits skal have en e-value
på 1e-10 eller bedre (mindre) for at vi regner det for
pålideligt. Dette passer også fint med, at disse hits er
enige om at der er tale om en histon.
Som udgangspunkt er der derfor 30 gode hits. Man kan dog godt
argumentere for at man skal passe på hits der kun har en
forudsagt funktion ("PREDICTED") - dem har vi
fire af. Der er også tre hits der hedder helt præcist det
samme ("histone H2A
[Homo sapiens]"), om her kunne
men evt. gå ind på link'et til sekvensen og læse hvad
der står - her kan man se, at der faktisk er mindre forskelle i
sekvens, og at de er stammer fra forskellige gener (/gene="HIST3H2A", /gene="HIST1H2AA",
/gene="HIST2H2AB).
Men lad os arbejde videre med alle 30 gode hits i de næste
analyser.
- Fra hvor mange gener stammer disse?
Svar:
Trick'et er her af en del af hits'ne er fra isoformer af det sammen
protein (se næste spørgsmål). Et protein med to
isofomer stammer stadig kun fra en enkelt gen.
Lad os fx kigge på følgende hits:
ref|NP_613075.1|
H2A histone family, member Y isoform 1 [Homo sa
149 7e-37
ref|NP_004884.1|
H2A histone family, member Y isoform 2 [Homo...
149 7e-37
Hvis man
gå ind og læser sekvens-entry'et for disse hits kan man
faktisk direkte se, at de stammer fra sammen gen (/gene="H2AFY", /gene="H2AFY").
Vi skal altså slå hits sammen, der dækker over
isoformer kodet fra sammen gen. Der er 5 "V" isoformer, 3 "Y" isoformer
og 2 "J" isoformer, så det bliver til 30 - 4 -2 -1 = 23 gener.
- Er der nogen af disse gener der har iso-former (altså
forskellige protein-produkter fra samme gen, typisk pga. alternativ
splejsning)?
Svar: Se
ovenstående.
- Bemærk: Noget af
sekvensen er markeret i gråt - det er lav-komplekse områder (dvs.
at de indeholder brudstykker af ofte forekommende sekvens, som ikke er
så meget værd for BLAST algorithmen). Man kan vælge
hvordan BLAST skal tage sig af disse område på den
første side. Hvis man vælger "low complexity" under
"filter", bliver disse områder slet ikke taget med når de
to indledende hits skal findes (default
i almindelig BLASTP!).
- Prøv at kører søgning igen (i et nyt vindue
eller tab) og denne gang slå lav-kompleksitets filteret til.
- Hvordan ser det nu ud med længderne af de alignments man
får?
Svar: De
bliver kortere i de tilfælde, hvor det lav-komplekse område
ligger i starten eller i slutningen af hit'et.
- Afsluttende
bemærkning: Vi har nu med BLAST fundet en række
homologe gener (eller rettere protein-produkter). Hvis man
ønsker at grave dybere, vil det næste logiske trin
være at trække den komplette sekvens for kandidat
proteinerne ud af databasen (og evt. ignorer isoformerne i
første omgang), for derefter at analysere dem med enten parvis alignment (for at finde
specifikke ligheder/forskelle mellem HTA2 og kandidaten) eller multiple alignment (kommende emne)
for at kunne analysere kandidaternes indbyrdes slægskabsforhold.
BLAST kan altså også være en god måde at
begynde at bygge et datasæt af sekvenser ud fra en kendt sekvens
(enten fra en anden organisme eller simpelthen bare et andet gen fra
samme organisme).
|