I have a issue to find out, is there any way that one can determine a group delay from frequency response of the FIR filter. Let's say I have I type FIR filter with order of 51. Frequency response obtained from freqz MATLAB's function is presented below. I know group delay is negative derivative of linear (in this example) phase response, but how can one determine that from the picture above.
The group delay is the negative derivative of the phase response as the OP has stated, and specifically for the delay of one clock sample, the phase will go linearly negative to $2\pi$ radians as the frequency goes from 0 to $f_s$
From the picture we see the phase is going approximately 800 degrees at a frequency of $.17\pi$ rad/sample (where $2\pi$ rad/sample is the sampling rate). So this would be equivalent to
$800/360*(2\pi)$ rad /$(.17\pi)$ rad/sample = $26$ sample delay
For linear phase filters, the group delay is half the order which would be 25.5 in this case (we can't resolve that from the graphic alone). The OP wrote that the order of the filter was 51 which implies 52 taps.