cannin/biostars_qa huggingface dataset:
Question: We can use edgeR to normalized raw read counts with TMM method, but how can we normalized TPM with the same method?

Answer: I don't quite follow your question. TMM is a method for normalizing the library sizes rather than a method for normalizing read counts. As the edgeR User's Guide says (page 15):

> normalization in edgeR is model-based, and the original read counts are not themselves
transformed.

Which way around is your question? Do you have TPMs and want to compute TMM factors or do you have TMM factors and want to compute TPMs?

If you are asking the first question, then no, TMM factors can only be computed from the raw counts, not from quantities such as TPMs or CPMs from which the library sizes have already been divided out. If you already have TPMs from some software package, then normalization has almost certainly already been applied, so I would be very wary about trying to re-normalize them unless you really know what you're doing.

If you are asking the second question then, yes, TMM factors can in principle be used to compute TPMs. In edgeR, any downstream quantity that is computed from the library sizes will incorporate the TMM factors automatically, because the factors are considered part of the effective library sizes. TMM normalization factors will be applied automatically when you use

```
CPM <- cpm(dge)
```

or

```
RPKM <- rpkm(dge)
```

in edgeR to compute CPMs or RPKMs from a DGEList object. I don't necessarily recommend TPM values myself, but if you go on to compute TPMs by

```
TPM <- t( t(RPKM) / colSums(RPKM) ) * 1e6
```

then the TMM factors will naturally have been incorporated into the computation.
Question: Hi all,

Could I please know if there is quickest way to extract reads from a huge fastq file to another. I already tried the following.

grep -A3 '1:N:0:' ORGAN1.fastq >ORGAN1.cleaned.fastq

but grep takes too long. Any oneliners from you are very much appreciated.

Thank you

Prakki Rama.

Answer: A faster way is to do this:

LC_ALL=C fgrep -A3 '1:N:0:' ORGAN1.fastq >ORGAN1.cleaned.fastq

Here is an explanation of why this is so much faster.
Question: This is probably a dummy question, but I haven't found an explanation and I want to be sure. When paired-end sequencing is done, there are usually 2 fastq files generated, one with "left" mates and the other with "right" mates.

Assuming that mates in the first file are read from the forward strand, are then mates in the second file read from the reverse strand?

If so, what is exactly written in the second file? A reverse complement of the read sequence (the corresponding sequence on the forward strand) or the pure reverse-strand sequence? And in which direction?

Answer: Illumina paired-end sequencing is based on the idea that you have initial DNA fragments (longer than your actual read length) and you sequence both its ends. On the Illumina chip, both ends of each sequence are amplified prior to actual sequencing using [bridging][1]. This approach results in two reads per fragment, with the* first read in forward* orientation and the *second read in reverse-complement* orientation. Depending on the inital fragment size and read length, these fragment can either overlap or not

For example, with 100 bp reads:

```
------[f-180bp]----->
----[r1]----->
<----[r2]----

------------------[f-350bp]----------------->
----[r1]----->
<----[r2]----
```

Therefore, the first fastq file will contain all "r1" reads, the second file all "r2" reads.

[1]: http://www.illumina.com/documents/products/techspotlights/techspotlight_sequencing.pdf
Question: Hi there!

I'm blasting (blastx2.2.28) my assembly to the annotated refseq. Below - the command line I used:

blastx -query assembly.fa -out Blast -db BlastCustomDB -num_threads 8 -evalue 1e-6 -best_hit_score_edge 0.05 -best_hit_overhang 0.25 -max_hsps_per_subject 1 -max_target_seqs 1 -outfmt 6

I specified quantity of HighcoringSegmentPairs (1), but still BLAST output contains two or more HSP for some queries.

Here's the output example (it goes like query ID, subject ID, identity %, length of HSP, mismatches, gap opens ... E-value, bit score):

TR60-c0_g1_i1	EKV52009.1	37.36	174	57	7	...	2e-14	74.3
TR60-c0_g1_i1	EKV52009.1	64.37	87	30	1	...	6e-10	59.7
TR61-c0_g1_i1	EKV51611.1	49.41	253	121	2	...	4e-47	161
TR64-c0_g1_i1	EKV51289.1	60.00	175	62	1	...	9e-66	208

After examination with 'uniq' command it became clear that all duplicated results for a query each time were aligned to a *single gene in different position* (what is exactly HPSs).

There's the way just [to remove multiple HSP][1], but maybe it's possible to change some parameters in blast run to get the right single-hit-single-HSP output directly?

Thanks in advance!

[1]: https://www.biostars.org/p/90188/

Answer: In NCBI BLAST+ versions below 2.2.29, the -max_hsps_per_subject parameter does not work correctly. This parameter was further replaced with -max_hsps, which works fine.

Download the latest NCBI BLAST+ (from [NCBI FTP][1]) and you will be fine:

blastx -query assembly.fasta -db blastdb.fa -max_hsps 1 -outfmt 6

Output:

TR60-c0_g1_i1	EKV52009.1	37.356	174	57	7	809	288	273	394	1.32e-18	74.3

[1]: ftp://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/LATEST
Question: I have a file with multi fasta format sequence

>lyrata FEN1
MGIKGLTKLLADNAPSCMKEQKFESYFGRKIAVDASMSIYQFLIVVGRTGTEMLTNEAGE
VTSHLQGMFNRTIRLLEAGIKPVYVFDGKPPELKRQELAKRYSKRADATADLTGAIEAGN
KEDIEKYSKRTVKVTKQHNDDCKRLLRLMGVPVVEATSEAEAQCAALCKSGKVYGVASED
MDSLTFGAPKFLRHLMDPSSRKIPVMEFEVAKILEELQLTMDQFIDLCILSGCDYCDSIR
GIGGQTALKLIRQHGSIETILENLNKERYQIPEEWPYNEARKLFKEPDVITDEEQLDIKW
TSPDEEGIVQFLVNENGFNIDRVTKAIEKIKTAKNKSSQGRLESFFKPVANSSVPAKRKG
NLICLVVSEVIPQYTGSQGTKCLLRQYKPAISNKVSSLHAFGWNSSKNIVSIRLFSLSWG
FNLSEIPESTTKGAANKKTKGAGGRKKK

>Arabidopsis thaliana FEN1
MGIKGLTKLLADNAPSCMKEQKFESYFGRKIAVDASMSIYQFLIVVGRTGTEMLTNEAGE
VTSHLQGMFNRTIRLLEAGIKPVYVFDGKPPELKRQELAKRYSKRADATADLTGAIEAGN
KEDIEKYSKRTVKVTKQHNDDCKRLLRLMGVPVVEATSEAEAQCAALCKSGKVYGVASED
MDSLTFGAPKFLRHLMDPSSRKIPVMEFEVAKILEELQLTMDQFIDLCILSGCDYCDSIR
GIGGQTALKLIRQHGSIETILENLNKERYQIPEEWPYNEARKLFKEPDVITDEEQLDIKW
TSPDEEGIVQFLVNENGFNIDRVTKAIEKIKTAKNKSSQGRLESFFKPVANSSVPAKRKG
NLICLVASEVIPNILDHRVPNACSGSINQPSAIKSPLCMLLDGTPPITLLALGCSLYPGL
PTFRFTWLEEIPESTTKGAANKKTKGAGGRKKK

>Brachypodium distachyon - II FEN1
MGIKGLTKLLAEHAPRAAAQRRVEDYRGRVIAIDASLSIYQFLVVVGRKGTEVLTNEAEG
LTVDCYARFVFDGEPPDLKKRELAKRSLRRDDASEDLNRAIEVGDEDSIEKFSKRTVKIT
KKHNDDCKKLLRLMGVPVVEAPGEAEAQCASLCKNHKAYAVASEDMDSLTFGSLRFLRHI
TDLSFKRSPVTEFEVPKVLEELGLTMDQFIDLCILSGCDYCENIKGIGGQRALKLIRQHG
CIEEVVQNLNNRFTVPEDWPYQEVRTLFKEPNVCTEIPDFQWTSVDKEGIVNFLAIENSF
SSDRVEKAVQKIKAARDRYSPGRVKLLTPVANLSGSITKKEPECVLGSSGQGLKSWSALQ
VCRSSSSDFRYRSYYSSKPLVLGRQSGFHGIPHAFSLI

and i want to extract sequence of defined region. these parameters are present in other file

lyrata 1 108
lyrata 147 234
Arabidopsis 1 108
Arabidopsis 147 234
Brachypodium 1 93
Brachypodium 132 219

I have asked similar question before and got an awk solution but it gives error. Please help me. i am a window user.

Answer: You can also use the [samtools][1] package. Samtools allows you to index a fasta file. Here is an excerpt from the [manual][2]:

> samtools faidx ref.fasta [region1[...]] : Index reference sequence in the FASTA format or extract subsequence from indexed reference sequence. If no region is specified, faidx will index the file and create .fai on the disk. If regions are specified, the subsequences will be retrieved and printed to stdout in the FASTA format.

Once you indexed your reference, you can query it.

For example:

samtools faidx reference.fasta
samtools faidx reference.fasta lyrata:1-108

Samtools is designed to handle nucleotides, so I am not sure that it will work with protein sequences.

To do it (in perl) for several positions (given in a file formated like in your question):

open(POSITIONS,"positions.txt");
while(){
chomp;
my ($seqName,$begin,$end) = split(/\s/);
open(SAMTOOLS,"/path/to/samtools faidx /path/to/reference.fasta $seqName:$begin-$end |");
while(my $line = ){
print $line;
}
close(SAMTOOLS);
}
close(POSITIONS);

[1]: http://samtools.sourceforge.net/
[2]: http://samtools.sourceforge.net/samtools.shtml
Question: I have the impression that FPKM as calculated by cufflinks doesn't correctly normalize for library size in case of unequal mapping rates between samples (e.g. in case of contamination). I get this impression from real data.

As a hypothetical example let's say that we have two libraries A,B, B 10 times the number of reads as A. One might expect, (I know it's a rough estimate) for a large number of genes in B to have 10 times the number of reads as in A, but mostly equal FPKM except for lowly expressed transcripts. Now assume A has ~100% mapping rate, and B only 50% (due to contamination, 40% will align to other references if known). Then, **many transcripts in B still seem to have 10X as many reads mapped as in A (maybe because number of mapped reads is not influenced linearly by contamination?)**, but the effective aligned library size that cufflinks sees is only 50% of the true library size. That would lead to an average 2-fold difference in median of all FPKM values.

Reference size is not taken into account for FPKM, is it? If I knew all contaminant genomes beforehand, I could just add them to the reference, and achieve a much more homogeneous mapping rate and better normalization between samples. On the other hand this might not be a fair comparison, as the contaminants would never have a chance of getting a FPKM assigned and will not influence the calculated average FPKM.

Is this correct from your experience, how to best compensate for it:

1. Add all possible contaminants to reference before alignment?
2. Use the sequenced library size instead of aligned library size?
3. Assume equal mapping rate for all samples (e.g. median, would be similar to median scaling of FPKM)?
4. Apply quantile normalization to FPKM data?
5. FPKM is already normalized, don't apply double normalization....

Answer: Well, some (including one of the original authors of the paper that introduced it) say that the RPKM/FPKM concept is flawed altogether.

Honestly I am not quite sure what to make of that statement. If indeed true then it is not just a casual observation. It would mean the ALL papers that use RPKM to compare gene expression are wrong. But then the observation was made in 2012, RPKM/FPKM are still heavily used. So there.

http://blog.nextgenetics.net/?e=51

and here: https://www.youtube.com/watch?v=5NiFibnbE8o
Question: I've wrote this multi-facted question to solve a specific coding problem, but also to bring together information that is covered partially across multiple existing biostars posts that others may find useful.

Part 1 - Edit the bam file
------

I have several bam files that have reads mapped to chromosomes that I'm not interested in i.e.

- 'chrUn_gl0000230' - over 40 of these
- 'chr9_gl000198_random' - over ten of these
- 'chr6_cox_hap2' - around 10 of these
- 'chrM'

The unique identifiers for the first three groups are 'chrUn', 'random', 'hap'.

I tried this to get rid of everything except **chr 1-22** and **chr X-Y** by cobbling together code from [here][1] and elsewhere:

samtools idxstats my.srtd.bam | grep -v chrUn* | grep -v *random | grep -v chrM | grep -v *hap* | xargs samtools view -b my.bam > my.only1to22andXY.srtd.bam

However, this eliminated only the **chrM** reads and the **chrUn** reads. I also tried this which gave the same result:

samtools idxstats my.srtd.bam | grep -v chrUn_* | grep -v chr*_random | grep -v chrM | grep -v chr*_hap* | xargs samtools view -b my.srtd.bam > my.only1to22andXY.srtd.bam

I suspect I'm having an issue with the * character but not sure why some chromosome reads are omitted and others are not.

What I really need is a one liner that **retains** only reads mapped **chr1-22** (I don't need the chr X-Y either). Is this possible? It doesn't matter if it involves removing what I don't need, or keeping what I do need. Ideally this would be in **bash**.

Part 2 - Edit the Header
------------------------

This leads to the next question of how to update the headers of the bam files. I need to omit the extraneous chromosome reads for a downstream application (ALTRE) as it crashes when they are present. There is a chance the package may also crash if the extra chromosome names are quoted in the header (even if the reads have been removed) and they read 0. Is it possible to clean the header to remove lines pertaining to extraneous chromosomes? Such as:

@SQ	SN:chr11_gl000202_random	LN:40103

Part 3 - Edit the fasta file
----------------------------

I'm aligning to hg19. I just occurred to me that the best course of action (particularly for future analyses) is to remove these extraneous chromosomes at source, i.e. from the fasta files or the genome index files that I use for the alignment. Is that a valid thing to do? Could anyone suggest how to go about this?

[1]: https://www.biostars.org/p/46327/

Answer: create a BED for chr1 to chr22 from the fasta index of your reference genome.

awk '/^chr[0-9]*\t/ {printf("%s\t0\t%s\n",$1,$2);}' human_g1k_v37_prefix.fasta.fai > your.bed

use this bed file to filter with samtools.

samtools view -L your.bed -o out.bam input.bam

---

another idea:

samtools view -o out.bam input.bam `seq 1 22 | sed 's/^/chr/'`
Question: Hi,

I am trying to search for the presence of couple sequences (around 400) each with a size of 23 bps,in different fastq files, while allowing 1-2 mismatches at maximum.
I am not sure if turning the fastq to a genome(transcriptome) would be a nice approach? I have tried making the fastq -> fasta -> building blast database -> running blastn, however it did not run as my query is not only one sequence.

> Example part of my query.file :

> ATTTTTCTGAAAAACCCCCTACGA

> AACAGGAAGTCAAAAAAAGCCAA

> AGGATTTTTTTTTTTCTGGGGACA

The output I am aiming to have is, for each read in my query.file, which of these sequences are having 100% (or having 1-2 mismatches) match in fastq file, and possibly where in the fastq file.

I would appreciate your suggestions! Thank you!

Answer: You can grab the fastq sequences containing these 23-mers with BBDuk like this:

bbduk.sh in=file.fastq outm=matched.fastq ref=23mers.fa k=23 hdist=2

"hdist=2" allows 2 mismatches; you can alteratively set that to 1 or 0. This does not tell you where the match is, but you can do that like this:

bbduk.sh in=matched.fastq out=masked.fastq ref=23mers.fa k=23 hdist=2 kmask=lc

That will convert the matched regions to lowercase.
Question: My RNA-Seq data is in format of fastq(ungzipped from fastq.gz format).I used STAR 2.5.3a mapping the reads with already indexed reference genome. It seems good. But I found the size of generated SAM file is strange. My original input fastq data is like 1.3-1.5 GB, but the SAM file ranges from 3.8 GB to 4.5 GB. Is that normal? If not, what is something wrong there?

Answer: A SAM file will typically be larger than a FASTQ file because, in general, it contains all the information of the FASTQ plus a lot of other information.

In addition, each FASTQ record may produce more than one alignment, hence you can see how it could easily grow to be much larger than the original data.
Question: Hi,

I have a large file of each chromosome for sequence data and I want to get the average score in 1kb bins?
Is there an easy way to do to this?

chr1 123 127 1890
chr1 208 210 90
chr1 215 220 50
chr1 225 230 100
chr1 250 300 800

I have tried bedtools, it will create the windows but then how do I get the average of the fourth value in that window? It also seems to be quite slow.

bedtools makewindows -w 1000 -g chromosome.txt > windows.txt

Any help would be much appreciated.

Answer: for fun, using awk+sqlite3, window size=100:

rm -f tmp.db && awk 'BEGIN {printf("create table T(C TEXT,S INT,E INT,V INT); BEGIN TRANSACTION;\n");}{printf("INSERT INTO T(C,S,E,V) VALUES(\"%s\",%s,%d,%s);\n",$1,$2,$3,$4);} END {printf("COMMIT; SELECT C,(S/100)*100 as G,((S/100)+1)*100,AVG(V) FROM T GROUP BY C,G;");}' input.bed | sqlite3 -separator $'\t' tmp.db && rm -f tmp.db

chr1	100	200	1890.0
chr1	200	300	260.0
Question: **How do you classify the miRNAs? What is the best way to classify miRNAs hypothetical using Blastn-short or RFAM? [doubts]**

Hello everyone. Let's me explain my doubts. I mean about 200 sequences which was predicted by mir-Bag (pre-microRNA) before. Then I used this pre-micros with the miRdup in order to obtain the microRNAs. Now, I want to know which have been listed above. Firstly, I did a `blastn -task blastn-short -evalue 10` (v. 2.2.30) against the mature miRNAs from mirBase.org. I wonder if you use this strategy and what options you usually use like:

```
-evalue 10
-penalty -4?
-word_size 4 or 7?
-ungapped?
-reward 1?
```

Because it is not neccesary to have a full coverage of the query, therefore there are some isomeRs with slightly mismatches. But if the seed region (2nt-7nt) is covered fully then it is the same family miRNA. Isn't it? Or I have to do the RFAM. Regards.

I attach one example and some options of blastn (v. 2.2.30.)

```
query= ENTRY327687_trypsinogen_354_553

Length=22
Score E
Sequences producing significant alignments: (Bits) Value

dme-miR-4941-5p MIMAT0020148 Drosophila melanogaster miR-4941-5p 24.3 0.21
dvi-miR-9538-3p MIMAT0035619 Drosophila virilis miR-9538-3p 22.3 0.83
mmu-miR-7677-5p MIMAT0029868 Mus musculus miR-7677-5p 22.3 0.83
mmu-miR-6922-3p MIMAT0027745 Mus musculus miR-6922-3p 22.3 0.83
gga-miR-6708-5p MIMAT0025818 Gallus gallus miR-6708-5p 22.3 0.83
mtr-miR2590j MIMAT0021334 Medicago truncatula miR2590j 20.3 3.3
mtr-miR2590i MIMAT0021333 Medicago truncatula miR2590i 20.3 3.3
mtr-miR2590h MIMAT0021332 Medicago truncatula miR2590h 20.3 3.3
tca-miR-3813-5p MIMAT0018648 Tribolium castaneum miR-3813-5p 20.3 3.3
mtr-miR2598 MIMAT0013301 Medicago truncatula miR2598 20.3 3.3
oan-miR-181c-3p MIMAT0007059 Ornithorhynchus anatinus miR-181c-3p 20.3 3.3

Query_29 2 TGGCGGTGAGCAGAATAATTG 22
20244 19 GGCGGTGAGCAG 8
34141 12 CAGAATAATTG 22
28818 4 GGTGAGCAGAA 14
26045 13 GGTGAGCAGAA 3
25597 11 GCGGTGAGCAG 21
21071 12 GCAGAATAAT 21
21070 12 GCAGAATAAT 21
21069 12 GCAGAATAAT 21
18383 10 GAGCAGAATA 1
13217 19 GCAGAATAAT 10
7532 10 TGGCGGTGAG 1

Lambda K H
1.37 0.711 1.31

Gapped
Lambda K H
1.37 0.711 1.31

Effective search space used: 4258254
```

Answer: I think BLAST is not the right tool, since it requires a word to be in exact match, you should try simple global alignment algorithm like needle of EMBOSS.

Rfam would be a good choice as well (or even better) but maybe you'll have a better success with the precursors, rather than the mature miRNAs since they contain more information.

Clustering of miRNAs is usually done using their seed sequences but this is not always true, some miRNAs can be in the same family (same seed) but have different targets since their 3' end changes specificity or vice-versa (different families - same targets).
Question: I am newbie to linux and NGS. Whle trying to build an index with hisat2 for my project, im getting following error:

"Error: could not open 2.fa
Total time for call to driver() for forward index: 00:00:00
Error: Encountered internal HISAT2 exception (#1)
Command: hisat2-build --wrapper basic-0"

the command and output of the terminal is as follows:

./hisat2-build /data/memona/reference/1.fa,2.fa,3.fa,4.fa,5.fa,6.fa,7.fa,8.fa,9.fa,10.fa,11.fa,12.fa,13.fa hisat_index

Settings:
Output files: "hisat_index.*.ht2"
Line rate: 6 (line is 64 bytes)
Lines per side: 1 (side is 64 bytes)
Offset rate: 4 (one in 16)
FTable chars: 10
Strings: unpacked
Local offset rate: 3 (one in 8)
Local fTable chars: 6
Local sequence length: 57344
Local sequence overlap between two consecutive indexes: 1024
Endianness: little
Actual local endianness: little
Sanity checking: disabled
Assertions: disabled
Random seed: 0
Sizeofs: void*:8, int:4, long:8, size_t:8
Input files DNA, FASTA:
/data/memona/reference/1.fa
2.fa
3.fa
4.fa
5.fa
6.fa
7.fa
8.fa
9.fa
10.fa
11.fa
12.fa
13.fa
Error: could not open 2.fa
Total time for call to driver() for forward index: 00:00:00
Error: Encountered internal HISAT2 exception (#1)
Command: hisat2-build --wrapper basic-0 /data/memona/reference/1.fa,2.fa,3.fa,4.fa,5.fa,6.fa,7.fa,8.fa,9.fa,10.fa,11.fa,12.fa,13.fa hisat_index

the hisat2-build file have all the permissions

-rwxrwxrwx 1 memona memona 2807 Mar 12 2016 hisat2-build

need kind help to resolve the issue

Answer: Hi,

You're missing one thing with your parameter expansion in bash. Try the following commnd line:

./hisat2-build /data/memona/reference/{1.fa,2.fa,3.fa,4.fa,5.fa,6.fa,7.fa,8.fa,9.fa,10.fa,11.fa,12.fa,13.fa} hisat_index

It should be better, as I think all your fasta files are in the reference directory.
Have a look at this conversation in [Stack][1]

[1]: https://stackoverflow.com/questions/5349718/how-can-i-repeat-a-character-in-bash#5349842
Question: I am planning to use STAR to map am RNA-seq I have to a mouse genome. I have never done that before. According to the [manual][1], I need to give a path to fasta files. When I looked into [Ensembl][2] database, I found several fasta files for different chromosomes.

I would appreciate any opinion, advice, or link on how to choose fasta files for mouse.

Thank you!

[1]: https://github.com/alexdobin/STAR/blob/master/doc/STARmanual.pdf
[2]: ftp://ftp.ensembl.org/pub/release-100/fasta/mus_musculus/dna/

Answer: In the Ensembl link you attached, there are files that are named something-something-"assembly", for example Mus_musculus.GRCm38.dna.primary_assembly - that is what you need
Question: What happened with miRecords mirna database website? The webpage is not available and I used it's data as part of an analysis in my paper. Now the reviewer can't see the page and asks for an explanation.

Edit:

Sorry, I forgot to post the link. The miRecords was suppose to be available at http://miRecords.umn.edu/miRecords.

Answer: I have left a comment on the PubMed site, and invited the team to please comment. http://www.ncbi.nlm.nih.gov/pubmed/18996891
Question: Hello,

I am interested in extracting information about the length of specific Reactome pathways (preferably using Bioconductor).

That is, given a pathway ID, I would like to know how many links it contains. Also, if possible, I would like to know how many links point to my Entity of interest, and how many links are there after it.

I am pretty new to bioinformatics and I hope I could make myself clear.

Thanks!

Answer: You can use [graphite package](https://www.bioconductor.org/packages/release/bioc/html/graphite.html)

For example if you want to know how many nodes and edges are in the first reactome pathway:

#load package
library(graphite)
require(igraph)
#get all reactome pathways for human
reactomedata=pathways(species = 'hsapiens','reactome')
#information about the first pathway in the list
reactomedata[[1]]
# number of nodes
length(nodes(reactomedata[[1]]))
# number of edges
dim(graphite::edges(reactomedata[[1]]))[1]

You can also create an igraph to calculate more topological indexes:

# builds a graphNEL
g <- pathwayGraph(reactomedata[[1]])
# create an igraph from the graphNEL object
g2=igraph.from.graphNEL(g, name = TRUE, weight = TRUE, unlist.attrs = TRUE)
# calculate diameter
diameter(g2)
# calculate transitivity/clustering coefficient
transitivity(g2)
# calculate outdegree
degree(g2,mode = 'out')
# calculate indegree
degree(g2,mode = 'in')

Obviously the topological indexes will change greatly if you have a directed edges, if you want to remove self loops and so on.