try commenting out [[
[don.git] / posts / cytoband_information_in_ncbi.mdwn
1 [[!meta title="Finding out Cytobands/Idiograms for assemblies"]]
2
3 In many organisms it is common to use
4 [idiograms](https://secure.wikimedia.org/wikipedia/en/wiki/Cytogenetics#Advent_of_banding_techniques)
5 or cytobands which provide information on approximately where
6 something is located on a chromosome in reference to the chromosome's
7 larger structure, or when exact locations are not required.
8
9 Until recently, I didn't know where [NCBI](http://ncbi.nlm.nih.gov)
10 kept their idiogram annotations, which made my
11 [[mirror of dbsnp|genetics/dbsnp_mirror/]] (which I use to annotate my
12 whole genome analyses) slightly less useful than it could have been.
13 But, after a bit of searching of NCBI's ftp site, I was able to locate
14 the file in the new
15 [movie directory](ftp://ftp.ncbi.nlm.nih.gov/genomes/MapView/Homo_sapiens/objects/current/initial_release): 
16 `ideogram_9606_GCF_000001305.13_850_V1`.
17
18 Then, a quick bit of work with SQL, I have the following schema:
19
20     CREATE TABLE idiogram (
21            chr TEXT NOT NULL,
22            pq  TEXT NOT NULL,
23            idiogram TEXT NOT NULL,
24                    -- I think these are related to recombination rates, but I'm not sure
25            rstart INT NOT NULL,
26            rstop INT NOT NULL,
27            start INT NOT NULL,
28            stop INT NOT NULL,
29                    -- I believe this indicates whether the band is black or white
30            posneg TEXT NOT NULL
31            );
32            
33     CREATE UNIQUE INDEX ON idiogram(chr,pq,ideogram);
34     CREATE UNIQUE INDEX ON idiogram(chr,start);
35     CREATE UNIQUE INDEX ON idiogram(chr,stop);
36
37 and an additional bit of SQL in my SNP annotation perl script:
38
39     SELECT CONCAT(chr,pq,idiogram) AS idiogram
40       FROM idiogram
41       WHERE idiogram.chr = ? AND idiogram.start <= ? AND idiogram.stop < ? LIMIT 1;
42
43 and some code:
44
45     sub find_idiogram {
46         my %param = @_;
47     
48         my %info;
49         my $rv = $param{sth}->execute($param{chr},$param{pos},$param{pos}) //
50         die "Unable to execute statement properly: ".$param{dbh}->errstr;
51         my ($idiogram) = map {ref $_ ?@{$_}:()} map {ref $_ ?@{$_}:()} $param{sth}->fetchall_arrayref([0]);
52         if ($param{sth}->err) {
53         print STDERR $param{sth}->errstr;
54         $param{sth}->finish;
55         return 'NA';
56         }
57         $param{sth}->finish;
58         return $idiogram // 'NA';
59     }
60
61 and viola:
62
63 | id         | chr |      pos | idiogram | ref | alt | orig_id    | gene   | [...] |
64 |------------|-----|----------|----------|-----|-----|------------|--------|-------|
65 | rs10000010 |   4 | 21618674 | 4p16.3   | T   | C   | rs10000010 | KCNIP4 | [...] |
66
67 idiograms for every SNP.
68
69
70 [[!tag genetics snp biology tech]]