The basic syntax of NeighborSearch is the following. If you tell me what you want to do exactly I can give you a more precise hint. Also, see the Bio.PDB FAQ (pdf).

The script basically loads the local pdb structure file molecule.pdb, creates a list of all atoms, gets the coordinates of the first atom, and prints a list of all resides within 5 angstrom of it (source).

from Bio.PDB import *
import sys
structure = PDBParser().get_structure('X', 'molecule.pdb') # load your molecule
atom_list = Selection.unfold_entities(structure, 'A') # A for atoms
ns = NeighborSearch(atom_list)
center = atom_list[0].get_coord()
neighbors = ns.search(center, 5.0) # 5.0 for distance in angstrom
residue_list = Selection.unfold_entities(neighbor_list, 'R') # R for residues
print residue_list

If you are interested in the neighbors of just a few atoms/residues, you might want to consider PyMOL since it has a nice set of GUI features that are more easily accessible.

here is a script that calculates the fnat using modeller - and it is a bit more straightforward to tease apart the ligand and receptor (it should all work automatically, no messy renumbering stages). it is however quite slow, but maybe you will find it useful ... !

from modeller import *
from modeller.scripts import complete_pdb
import sys,os
import math
import numpy as np

def calculate_fnat(model,cut_off):
    fnat = {}
    for a in model.chains:
        for b in model.chains:
            if a != b:
                sel_a = selection(a)
                sel_b = selection(b)
                for ca_a in sel_a:
                    ax = ca_a.x
                    ay = ca_a.y
                    az = ca_a.z
                    for ca_b in sel_b:
                        bx = ca_b.x
                        by = ca_b.y
                        bz = ca_b.z

                        dist = math.sqrt(((ax-bx)**2)+((ay-by)**2)+((az-bz)**2))

                        if dist <= cut_off:
                            a_res  =  ca_a.residue
                            b_res  =  ca_b.residue
                            if (int(b_res._num)+1,b.name,int(a_res._num)+1, a.name) not in fnat:
                                if (int(a_res._num)+1,a.name,int(b_res._num)+1, b.name) not in fnat:
                                    fnat[int(a_res._num)+1,a.name,int(b_res._num)+1, b.name] = dist 
    return fnat

# initialise the environment .
env = environ()
env.libs.topology.read('${LIB}/top_heav.lib')
env.libs.parameters.read('${LIB}/par.lib')

#calculate native contacts
m1 = complete_pdb(env, 'native.pdb')
rmsd_sel = selection(m1).only_atom_types('CA')

fnat = calculate_fnat(m1,5)

filenames = []

for pdb_file in [l for l in os.listdir("./") if l.endswith(".pdb")]:
    print pdb_file
    filenames.append(pdb_file)

results = []    
# calculate model contacts                      
for fname in filenames:
    m2 = complete_pdb(env, fname)

    fnat_mod = calculate_fnat(m2,5)   

    count = float(0)

    #compare native contacts to model contacts
    for f in fnat:
        if f in fnat_mod:
            count += 1

    score = count/len(fnat)

    results.append([fname,score])

#output results
f1 = open('data_4.txt','w')
f1.write('fname\tfnat\n')    

for line in results:
    f1.write('%s\t%s\n' % (line[0],line[1]))

f1.close