<p>Here's a pure python solution to process your example. It should be straightforward to loop over all the rows of your vcf files</p>

<pre><code>my_line='16 32328771    esv3638416;esv3638417;esv3638418;esv3638419;esv3638420  A   &lt;CN0&gt;,&lt;CN2&gt;,&lt;CN3&gt;,&lt;CN4&gt;,&lt;CN5&gt;   100 PASS    AC=27,84,3220,20,26;AF=0.00539137,0.0167732,0.642971,0.00399361,0.00519169;AN=5008;CS=DUP_gs;END=32370968;NS=2504;SVTYPE=CNV;DP=26381;EAS_AF=0.0099,0.0188,0.5427,0.004,0.0089;AMR_AF=0.0101,0.0058,0.4063,0,0.0014;AFR_AF=0.0015,0.0371,0.8374,0.003,0.003;EUR_AF=0.003,0.0089,0.5189,0.005,0.005;SAS_AF=0.0051,0.0031,0.7791,0.0072,0.0072;VT=SV  GT  3|3 3|1 0|3 3|3 3|3 3|0 0|3 3|3 0|0 3|3 1|3 '
my_line=my_line.strip() # remove trailing whitespace
my_fields=my_line.split()
names=my_fields[4].split(',') # extract CNV names
names=[x[1:-1] for x in names] # remove '&lt;&gt;' symbols
print('\t'.join(names))
for sample in my_fields[9:]:
    genotypes=[0]*(len(names)+1) # initialize array of zero counts, including reference  
    for allele in sample.split('|'):
        genotypes[int(allele)]+=1
    print('\t'.join([str(x) for x in genotypes[1:]])) # print all columns except reference
</code></pre>

CC BY-SA 3.0
2018-01-04T20:45:47.880Z
Thanks @heathobrien, this may work but in the example above I just showed genotypes for a few individuals, I have nearly 900 individuals and It is not convenient to copy all genotypes in your code in "my_line" part! Do you have any alternative to be able to fix this?
Anna1364
CC BY-SA 3.0
2018-01-04T21:23:44.640Z
You'll probably want to include the code I wrote as part of a larger script that reads your genotypes file and writes the output to another file rather than printing it to the screen. This isn't a script writing service, so take a look at an python I/O [tutorial](https://www.tutorialspoint.com/python/python_files_io.htm) or the documentation for the csv module (or pysam)
heathobrien
CC BY-SA 3.0
2018-01-03T10:28:41.517Z
2018-01-03T10:28:41.517Z
heathobrien
<p>I have a bunch of files (vcf format) each containing site and sample annotations of a copy number variant (following example is part of one of my files, I am showing 11 individuals but they are nearly 900 individuals)</p>

<pre><code>16  32328771    esv3638416;esv3638417;esv3638418;esv3638419;esv3638420  A   &lt;CN0&gt;,&lt;CN2&gt;,&lt;CN3&gt;,&lt;CN4&gt;,&lt;CN5&gt;   100 PASS    AC=27,84,3220,20,26;AF=0.00539137,0.0167732,0.642971,0.00399361,0.00519169;AN=5008;CS=DUP_gs;END=32370968;NS=2504;SVTYPE=CNV;DP=26381;EAS_AF=0.0099,0.0188,0.5427,0.004,0.0089;AMR_AF=0.0101,0.0058,0.4063,0,0.0014;AFR_AF=0.0015,0.0371,0.8374,0.003,0.003;EUR_AF=0.003,0.0089,0.5189,0.005,0.005;SAS_AF=0.0051,0.0031,0.7791,0.0072,0.0072;VT=SV  GT  3|3 3|1 0|3 3|3 3|3 3|0 0|3 3|3 0|0 3|3 1|3     
</code></pre>

<p>I want to count each allele of CNV in each individual(sample) if there are 0, 1 or 2 copies!
this is my desired output (rows are individuals):</p>

<pre><code>CN0 CN2 CN3 CN4 CN5
0   0   2   0   0
1   0   1   0   0
0   0   1   0   0
0   0   2   0   0
0   0   2   0   0
0   0   1   0   0
0   0   1   0   0
0   0   2   0   0
0   0   0   0   0
0   0   2   0   0
1   0   1   0   0
</code></pre>

<p>Any idea how to do that?</p>

CC BY-SA 3.0
2018-01-02T10:04:04.683Z
Cold you explain your file format a bit? How do we get from the VCF line you show to the table of your desired output?
terdon
CC BY-SA 3.0
2018-01-03T00:20:44.670Z
just edited my post to clarify that it is a vcf form file. I used vcftools --poition option to extract the position that I needed! and how to create the table is exactly the question that I have .....!
Anna1364
<python><perl><cnv>
how to count the number of each allele of CNVs?
# how to count the number of each allele of CNVs?

(Asked by: username_0 on 2018-01-01)

I have a bunch of files (vcf format) each containing site and sample annotations of a copy number variant (following example is part of one of my files, I am showing 11 individuals but they are nearly 900 individuals)
16  32328771    esv3638416;esv3638417;esv3638418;esv3638419;esv3638420  A   ,,,,   100 PASS    AC=27,84,3220,20,26;AF=0.00539137,0.0167732,0.642971,0.00399361,0.00519169;AN=5008;CS=DUP_gs;END=32370968;NS=2504;SVTYPE=CNV;DP=26381;EAS_AF=0.0099,0.0188,0.5427,0.004,0.0089;AMR_AF=0.0101,0.0058,0.4063,0,0.0014;AFR_AF=0.0015,0.0371,0.8374,0.003,0.003;EUR_AF=0.003,0.0089,0.5189,0.005,0.005;SAS_AF=0.0051,0.0031,0.7791,0.0072,0.0072;VT=SV  GT  3|3 3|1 0|3 3|3 3|3 3|0 0|3 3|3 0|0 3|3 1|3     

I want to count each allele of CNV in each individual(sample) if there are 0, 1 or 2 copies!
this is my desired output (rows are individuals):
CN0 CN2 CN3 CN4 CN5
0   0   2   0   0
1   0   1   0   0
0   0   1   0   0
0   0   2   0   0
0   0   2   0   0
0   0   1   0   0
0   0   1   0   0
0   0   2   0   0
0   0   0   0   0
0   0   2   0   0
1   0   1   0   0

Any idea how to do that?

**Tags:** python, perl, cnv
> **Comments:**
> **username_1**: Cold you explain your file format a bit? How do we get from the VCF line you show to the table of your desired output?  
> **username_0**: just edited my post to clarify that it is a vcf form file. I used vcftools --poition option to extract the position that I needed! and how to create the table is exactly the question that I have .....!  

---

## 1 answers

### Answer by username_2 on 2018-01-03. Score: 1 | Accepted answer

Here's a pure python solution to process your example. It should be straightforward to loop over all the rows of your vcf files
my_line='16 32328771    esv3638416;esv3638417;esv3638418;esv3638419;esv3638420  A   ,,,,   100 PASS    AC=27,84,3220,20,26;AF=0.00539137,0.0167732,0.642971,0.00399361,0.00519169;AN=5008;CS=DUP_gs;END=32370968;NS=2504;SVTYPE=CNV;DP=26381;EAS_AF=0.0099,0.0188,0.5427,0.004,0.0089;AMR_AF=0.0101,0.0058,0.4063,0,0.0014;AFR_AF=0.0015,0.0371,0.8374,0.003,0.003;EUR_AF=0.003,0.0089,0.5189,0.005,0.005;SAS_AF=0.0051,0.0031,0.7791,0.0072,0.0072;VT=SV  GT  3|3 3|1 0|3 3|3 3|3 3|0 0|3 3|3 0|0 3|3 1|3 '
my_line=my_line.strip() # remove trailing whitespace
my_fields=my_line.split()
names=my_fields[4].split(',') # extract CNV names
names=[x[1:-1] for x in names] # remove '<>' symbols
print('\t'.join(names))
for sample in my_fields[9:]:
    genotypes=[0]*(len(names)+1) # initialize array of zero counts, including reference  
    for allele in sample.split('|'):
        genotypes[int(allele)]+=1
    print('\t'.join([str(x) for x in genotypes[1:]])) # print all columns except reference

> **Comments:**
> **username_0**: Thanks @heathobrien, this may work but in the example above I just showed genotypes for a few individuals, I have nearly 900 individuals and It is not convenient to copy all genotypes in your code in "my_line" part! Do you have any alternative to be able to fix this?  
> **username_2**: You'll probably want to include the code I wrote as part of a larger script that reads your genotypes file and writes the output to another file rather than printing it to the screen. This isn't a script writing service, so take a look at an python I/O [tutorial](https://www.tutorialspoint.com/python/python_files_io.htm) or the documentation for the csv module (or pysam)  

---

